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This thesis presents an analysis of transonic potential 
flow through an oseillating unstaggered thin plate cascade. 
ιν Sol ocacion technique is used involving the superposition 
of adjacent blade "isolated foil" potentials with inter- 
ference potentials of unknown strength. Imposition of 
flow tangency requirements leads to integral equations for 
the unknown source distributions of the interference 
potentials. Results presented include the interference 
coefficients for the steady and oscillating cases. 

The steady case is extended to the determination of the 
potentials along the blade and compares favorably with a 
parallel solution using a Laplace transform approach to 


the sonic wind Tunnel wall interference problem. 
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I. INTRODUCTION 


Considerable advances are still being made in the 
development of aircraft turbopropulsion systems. Іп the 
Em ropr higher thrust to weight ratios, one must look to 
higher stage loading as well as significantly increased 
turbine inlet temperatures. From the compressor and fan 
Houenoepoint of view, this translates primarily into increased 
Mis peeds. Moving toward better blade design, the inves- 
raeation of transonic and supersoniec relative flows is 
becoming more important not only for increased performance 
ИШК а СО Гог fiutter prediction. Ав shown by Pratt and 
Whitney, various flutter problems have already been encoun- 
tered [Fig. l as extracted from Ref. 1] and more such 
problems can be expected in the future. 

Recognizing that the actual three dimensional flow 
fields encountered in a transonic stage are complicated by 
Mime tects Of rotation, finite blade thicknesses and cambers, 
boundary layers, and spanwise Mach distributions, it is 
excessively difficult to model these effects simultaneously. 
One can hope that eventually a complete analysis will 
become practical, but at this time the problems must be 
paced piecemeal. One simplified technique historically 
used in modeling such flows is the two dimensional cascade. 

The following investigation was concerned with transonic 


cascade flow. This problem had previously been considered 


E 
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by Hamamoto [Ref. 2] and Gorelov [Ref. 3]. Also of impor- 
tance is the related problem of transonic wind tunnel 
interference. Here the studies of Sandeman [Ref. 4] and 
Murasaki [Ref. 5] must be mentioned. For a more general 
ШУ КООП ОГ this problem, reference is made to Berndt 
[Ref. 6] and Goethert [Ref. 7]. 

Hamamoto based his analysis on the linearized oscilla- 
tory transonic flow equations.  Gorelov attempted to include 
the effect of blade thickness by approximating the nonlinear 
transonic perturbation equation. This approximation, 
however, is valid only for low supersonic flows. In the present 
Шімессірасіоіп, а transonic flow approximation was intro- 
duced based on previous work by Oswatitch [Ref. 8] and Teipel 
[Ref. 9]. The technique permitted a general formulation 
for both the steady choke flow problem (in conjunction with 
Spreiter's local linearization technique [Ref. 10]) and 
the oscillatory problem. The resulting boundary value 
problem was solved in a manner similar to that used by 
сео Line proolem formulation, solution technique, 

Me rical evaluation of the recurring integrals, and 
Pepresentative results are described in the following 


Sections. ~ 
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ІІ.  UNSTEADY TRANSONIC FLOW THEORY 


A. DEVELOPMENT OF THE GOVERNING TRANSONIC FLOW EQUATIONS 


Ue 





--Ше- 
0 2C A 


Figure 2.  COORDINATE SYSTEM 


η το lena rr L l Ow of an inviscid adiabatic perfect 
gas along the positive X axis which is perturbed by a thin 
body (or thin bodies) executing small transverse unsteady 


mocions. ithe perfect gas equation is applicable 
p = pRT (SO 


where: 
p = pressure 
p = density 
T = absolute temperature 


R = gas constant 
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вгпегпоге the conservation conditions for mass, momentum, 
 епетру тизс hold. Denoting the velocity components in 
the X and Y direction respectively as U and V, the conserva- 


tion of mass may be expressed as 





e + 90р) m ӘУЕ) = (axo) 


Conservation of momentum is expressed 


ЕШ aD δρ 
Do "ay БІ ЗЕ 
(II-3) 
DV. 39V, gàV ,. yv. 1 a3P 
Dt 3t ' o T Voy o ЭҮ 
Conservation of energy is 
DA EE W^, o M Ue ορ (II-H) 
Dt 2 pe `X oY 
ОР 
Ds . 
5-5 0 (II-5) 
where 
SE = substantial derivative 
e - energy for the individual particle 
S = entropy for the individual particle 
TEES 


T5 





ο ο аса or course, exelude finite discontinuities 
(shoeks) in the flow field. 

систе иот ап important role is played by certain 
functions which reduce the number of dependent variables. 
Nuscoften the velocity potential is introduced. Its 
ек сепсе depends on the eondition of irrotationality of 
τ τίου. For transonic flows over sufficiently thin 
bodies (small perturbation assumption!) it is usually 
possible to work with the velocity potential since any 
ПЕША су introduced by occurring shocks is of at least 
Unara order in the Е tion velocities. It follows 


Esc the velocity can be expressed as 


90 
U = — = 6 
` (II-6) 
29%. 
V AY Фу 


i Do - 2 


af = (с^; (Ir-8) 
s=constant 
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and that 


DOM MED. D 
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FU 


uu 


сї 


(у 


Therefore the eontinuity equation (II-7) becomes 


1 
D RD ES 


А ΥΕ ο {1510} 


Integration of the Euler equations along a streamiine gives 


the unsteady Bernoulli equation 


2 
W αμ. 


where F(t) = 502 luco unmMhorm parallel flow at infinity. 


Taking the partial of (II-11) with respect to t, one has 


2 
1 "q: W = 
o Pp = 560% + 5) Ми) 
juserting (11-3) апа (11-12) ѓог Р, Р, апа Р, іп (11-10) 


ХУ E 
and utilizing (II-6) one obtains the following partial 


differential equation in 0. 


° ° 6 < 26 à 
(1----)Ф + e A car o 
42 XX 4° уу ac ху 
20. 24 1 
= ay Ole = ae Sn το (11-13) 
а a a 
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Ihe usual way to simplify this equation is to neglect 


БЕ product terms obtaining 


= M 
(AMI un 2 L = LM - 0 
УУ со а une qe 


which is the well known unsteady linearized potential 
equation valid for the supersonic and subsonic regimes. 
1ου the case near M = 1, such a "naive" lineari- 
Шеш саппос, in general, be carried out. Consider the 
following analysis for the steady case. 


Define the perturbation potential, 4, such that 


Ó = UX + 9. (11-15) 


ἵ υπ steady case, the Bernoulli equation (II-11) is 


+ Е = constant . (а 07) 
Now since 
ГВ = 0 Ше гр, 
0 eM 
and 
а^ = (55 - YR T 


S 


Equation (II-16) becomes 
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- 


2 
+ τ - constant (II-17) 


ΙΙ, 


 сопссапг і5 еуа!цасей by noting that at the critical 


Mach number 


and (II-17) becomes 


x E 
s Y. oe Y= μ΄ EAD 


With (II-18) and (II-15) in (II-13), one has 
ср NS ү+1 2 
cae Ue) + фуу = ο ο... 5 ф, 
2 $ 1 + 26 [0 ф + фф 1] 
у I о V х'у 
= а +1 2 
ο ο ο πι ο 
ο. ο ne onedd pre yr ously; únder the Prandtl assump- 
tion, one would neglect the product terms on the right-hand 
ur το However. for M = 1, а поге сагегиі апа1уз15 


is in order. To this end one introduces the Prandti-Glauert 


transformation: 
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ф = єф(Х,Ү) 7220) 


where X 


= % 

Ү- Ү/1-м2 | aa) 
Us 

M = = Qu 


Here, one requires ọọ to be of order unity, therefore є 
assumes the order of magnitude of 6. 
Equation (II-18) becomes 
Do ο — — 
Eo (1-M )a. = £ (y+1)U 63595 


= 2. 2 = 
εϕ-- (1-M_)a, + E 


2 2 
tE a S 


“Ө- 
ее 


2*4. Q-MDU SL * 2e) Leer (IM) 


+ 


2 


++ 


Drs MAE LEM C - 
e^ (1-M^) (v Uy % к-15-(1-М) 7фууфұ 


2 
3 (Y-1) оъ 
MES mE Lu (II-23) 


> 


Now the terms on the RHS may be neglected for suitably 
A oyes trate the validity of such a measure, 
one may compare the largest term on the RHS with a term 
Mon tnem THS: This term is the first term of order P the 


remaining terms being of order c? ОГ Επ). Consider 
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the ratio 


ПЕСИ, %- 
Б = ———— ки 





EN cleapr that the nonlinear term becomes of higher 
order for the subsonic and supersonic cases but must be 
sped for traneonic flows. This consideration can be 
generalized to unsteady transonic flows thus leading to 
Piemeransonic small perturbation equation as the governing 
equation: 
> ма M, 
[ας (VFL) gr e a ὀχε “5 
тс, 
A ragorous justification for this equation was first given 
by Lin, Reissner, and Tsien [Ref. 11] and by Landahl [Ref. 12]. 
DEM berEssupiISCcation of this equation is possible if 
nero ransonic flow is sufficiently unsteady. Although this 
can be shown by mathematical order of иаа considera- 
ο ος το JS name references, it is instructive to consider 
the physical aspects of the transonic wave propagation 
problem. Figure 3 shows an acoustic source in transonic 
flow. It may be seen that the parts of the wave fronts 
moving upstream, referred to as receding waves, move away 
very slowly from the source. The downstream or advancing 
waves move with a velocity approximately equal to twice the 


speed of sound. 


2 








Figure 3. TRANSONIC WAVE PROPAGATION 


Bene now a continuous spatial distribution of distur- 
bance sources, it is readily understood from Fig. Y that the 
mecedine waves will tend to cancel provided the wavelength 
(given by 2m|a-U|/o , where ω is the circular frequency of 
the disturbance) is very short compared tobody reference 


πο ο ο However, if the sources oscillate very slowly, 





Source À and B with Respective Waves 


m a — 


The Resulting Wave 


Figure 4. WAVE CANCELLATION 
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pressure waves of the same sign will have sufficient time to 
η εγαος and possibly = ӘНШІСІ ро finite amplitudes. 
erzzreusiy, these receding waves may create large phase lags 
between motion and pressure on the airfoil causing severe 
Beroeclastic instability problems. 

νι. τπτ tne Conclusion, therefore, that sufficiently 
umsteady transonic flows сап be analyzed by linearized 
theory, Eqn. (II-14), whereas only mildly unsteady transonic 


flows must be investigated using Eqn. (II-25). The condition 


for linearization is generally stated as [Ref. 12] 


E SO ene 


where |1- M is the largest deviation of the Mach 


LOCAL! 
number from unity in the flow and K is the reduced frequency 


(K = wC/U_). 


B. DEVELOPMENT OF THE BOUNDARY CONDITIONS 

To completely describe the problem, proper boundary and 
initial conditions must be specified. The most important 
OMS Mare: 

i The flow on the surface of the body must be tangent 
ο any Jrnstarnt Of time. 

2. Sommerfeld's radiation condition must be satisfied 
A Sano ro time, 1.e., waves propagate away from 
sec csr qisturbanee. 

Further, certain edge and wake conditions must be 
satisfied, e.gp., smooth flow at trailing edges must be 


ес аиа иеа (Кисса condition). 
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The general flow tangeney condition can be stated as 


ο σης. If the equation of the surface of a body which is 


moving in a time dependent manner is given by 
Βλ χο) τρ (ІІ-26) 


ШІ Гле low tangeney condition requires that the total 


s (i ve Of the function F is zero, i.e., 
= = A = -- 
Е U В L V E 0 . MIS 


comes der an airfoil that is placed in a Cartesian coordinate 
man such a way that if the airfoil is restricted to 
EN ompiicude motions owt of the Y = O plane, then its 


шве апа lower surface can be described by 


DNE Y — Y (Xt) = 0 
u u 
(ІІ-28) 
E = Y — Y (X,t) = 0 
erore, from the boundary condition one has 
V = Yu, + Ὁ "e at Y = Y (X,t) 
(II-29) 
V = Үр ma Yr. at Y = Y (X,t) 
t X 
One now introduces the small perturbation assumption, 
Ó = UX + ¢ C 
where the disturbance velocity components are 
u = 0 
á (11-30) 
у = фу 
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auch that 


Б G" m 
V = v ΠΠ) 


11 
СЕ; 
a: 
= 


ШЕР U,V << UL 


In (II-29) the term, u Ү, > сап then be neglected compared 
me 
απ much larger term, Ue 00а оороо о оојоел о је = 


then simplify to 


μα у. πο νι να ο) 


(11-33) 
and v=Y ey Тор Ү - Yr (X,t). 


For small amplitude oscillations, (II-33) can be further 
simplified by expressing the flow tangency condition in the 


ре = 0. Taylor series expansion gives 


3v(X.0* ,t) 


л 
v(X,Y,,t) = νο ο 20) is Ya Doo o ος т» 


эм OS t) 


and | v(X,Y, ,t) — το ο 5) t Yr ovas ШІ PE 


Using the same argument aszabovesa che produet terms can 
be neglected so that the flow tangency condition can be 


stated as 


! 
O 


v = Y LEES ШЫГ Жү 
(11-24) 


II 
O 


v = Y + U Y TO CT 


For a thin foil performing a time dependent motion about 
the Y = 0 plane, the airfoil thickness terms can be dis- 


regarded by similar arguments. The above conditions become 


e5 





B 
Q2 
aJ 


y 55 for Y=0 (11-35) 


а? 
8 
а? 
Ps 


where h(X,t) describes the time dependent chord line of the 


blade. 


C. THE DERIVATION OF THE UNSTEADY EQUATIONS 


Define the following nondimensional terms. 








yu. X 
ἘΞ E 
Y = š 
(ІІ-36) 
pue 
f. 0 
ES U C 
qt OM. 
08 
— _ h 
h = š 


5ος the chord of the foil. Equations (II-25) 


and (II-35) may now be written in nondimensional form: 


2-- 


r 


(ІІ-37) 


ho 22 += T 2 
[1 - M. - Gau = a c M 


and ο ο ο. (ІІ-38) 


Since the main interest of this thesis centers оп the 
determination of the aeroelastic stability of a foil, only 


small perturbations from the steady state need be studied. 
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If one restricts oneself to such small amplitude perturba- 
τ... 15 possible to write the perturbation potential as 
ӨРШІП of a steady state and a time dependent potential. 
SOS in the casesof a harmonic oscillation, one may describe 
ri Lime dependent quantities as the product of a complex 


amplitude function and a harmonic time factor. Therefore 


S@emoetines the potentials 9, and y such that 


= iKt 


ФСХ,Ү,6) - %1(Х,Ү) + p(X,Y,K)e (11-39) 
Ehe foil chord function and vertical velocity are 
nowt) = h, (X) + К куе Gero) 
V(X,b) = v,(X) 8 v OX I)De ^ (11-41) 
оС 
where K = U (тада) 
It may be shown that (II-37) becomes 
[1-M2-m<(y+1)¢,—]¢,— + dy = 0 (11-13) 
о “co ex jo x x 1уу 


which is the steady state equation, and the unsteady equation 
is | 
2 2 9 2:2 
а с (9)-%-) + MK v 


- 


ΤΕ " 
| 2iK sum er ( 


The boundary equation (II-38) becomes 
v (X) = 210 Eo 


and v(X,K) = ho(X,k) + ikh(X,x) (II-46) 


e 





ΜΚ ОТ surface equation (11-43) must satisfy the 


boundary equation 


ο ο. τπτ est satisfy (IT-46). 


ШӘ сопуепзепсе, the overbar and hat will now be dropped 


and the steady and unsteady equations with boundary equa- 


tions become 


[1- 
subject to 

val 
and 

2 
ΙΙ. 

ο σου to 

V 


Memes that (II-l7) 


(II-49) is linear 


Мы - ME(V#L) Oy Ia + Фуу = 0 (II-17) 
ку at Y=0 uu ο 
- ме (үзі (9,29, + MÊK Y 

- 2iKMcy, + (uy = 0 (11-19) 
ο πο ο στο (11-50) 


is a nonlinear equation in 91. Equation 


in U but the existence of the variable 


coefficient 915 cons der abi ооо сарев 105 solution. 


ПО ΠΙΕΙ REVIEW OF POSSIBLE APPROACHES TO 
UNSTEADY TRANSONIC FLOW ANALYSIS 


The following is a presentation of various approaches to 


the solution of equation (II-49). 


1. Method of Teipel and Hosokawa 


An approximation of equation (II-49) was proposed 


by Teipel [Ref. 9] and Hosokawa [Ref. 13] for oscillating 
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τ: η 5οπῖς flow (M. = 1). For steady sonic flow, 
Oswatitch [Ref. 8] introduced a linearized approximation by 


setting 
(Ү%1)%4,, = À = constant, E NU +; (II-51) 
allowing (11-13) to be written in the easily solved form 


ф = x E = 0. (11-59) 


lyy 


Teipel and Hosokawa generalized this procedure to unsteady 
sonic flow by approximating equation (II-49) with (neglecting 


the term in ον) 


2 


-À um ΞΕ Voy - 21Ку, + K p= 0. (ІІ-53) 
m iming Q 
О- А 4 21К, τ 5η, 
amenas from (11-46) 
2 P 
Voy 7 Q Ye πο νο, (11-55) 


This equation allows the study of sonic flow past a 
single oscillating airfoil in unbounded flow over the entire 
frequency range. It contains the steady solution (Oswatitch 
approximation) as a special case. For large K (Refer to 
page 23) the term containing Pax becomes negligibly small; 
Gnereioreswretention of A at high reduced frequencies 
becomes unnecessary. Thus, when one takes À = O, the 


Teipel-Hosokawa theory also contains the linearized theory. 
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Local Linearization Sc NM Ue o 


Steady Transonic Flow 
It was shown that the Oswatitch approximation of the 


steady transonic small perturbation equation is 


ф - A oy = 0, A= CONS > 0 (ІІ-52) 


Eyy 


solution to this equation for sonic flow past a single 


БЕНШГОЛІ is 





2 
λΥ 
$, (X,Y) = - E f z : 179 ul (II-56) 
TA Χ-- 


МИГЕ ав interest is the velocity on the profile which 


is given by 


x 
EE а (EDU 
-1 VX-S 


ae 
PS 


$, 050) = - 2. 
Recalling the definition of А, (II-51), Spreiter argued that 
Er οριρίπαι variable, Φα» could be reintroduced into 
equation (11-57) thus linearizing the equation locally. 
This leads to an ordinary nonlinear differential equation in 
Py which can be solved analytically. The agreement with 
experiment is remarkably good [Ref. 10]. Its application 
to the choked wind tunnel wall interference problem was 
given by Sandeman [Ref. 4]. It is peered that an extension 
of this technique to the unsteady case may prove fruitful. 
3. WKBJ - Method 
Rubbert [Ref. 14] applied the .WKBJ method to this 


problem. Differentiating the axisymmetric counterpart of 
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equation (II-49) for M, = 1, multiplying by $4, and applying 


the Hankel transform, the equation ean be reduced to the 


form 


2 E 
ο ο ο D (11-58) 


where N is & large parameter. Unfortunately no inversion 
pue necresulrting transform integrals has been found. 
l. Method of Geometrical Acoustics 

A third approach has recently been proposed to obtain 
E Eenvekness effects for oscillating wings in transonic 
flow [Ref. 15]. This procedure begins with a given non- 
uniform steady transonic flow field which has been deter- 
mined either analytically or experimentally (In the latter 
case an analytical expression must be obtained which matches 
the measured pressure coefficient.). Superimposing upon 
ШЕЕ Low field small unsteady perturbations, the path and 
transmission time of these acoustic perturbations are 
Beccemened by the methods of geometrical acoustics. Again, 
only preliminary calculations have been carried out to date. 

5. Unsteady Transonic Flow with Shock Waves 

pic ci cr ensionle flows around an airfoil exhibit 
locally supersonic regions which are terminated by shock 
waves. Often, the interactions of small unsteady perturba- 
tions with these shock waves may be responsible for aero- 
dynamically induced instabilities. Eckhaus [Ref. 161 made 
the problem mathematically tractable by assuming a normal 


shock extending from the surface to infinity with uniform 


E 





ου ὃν supersonic flow in front and uniform steady subsonic 
IE behind the shock. ο sns a small perturbation 
analysis, i.e., superimposing small time dependent perturba- 
Bons upon this idealized Svecdy asocian LON. me was able 
memoovain equations for the interaction of the sound waves 
ШІ Ге shock. However, only a few numerical results were 
obtained. 


For the following study, the Teipel-Hosokawa 


approach was used. 


ge 





III. TRANSONIC UNSTEADY CASCADE FLOW PROBLEM 


A. PROBLEM FORMULATION FOR SONIC FLOW PAST 
OSCILLATING UNSTAGGERED CASCADES 


pa 


Uso Πο 
Eu 


Figure 5. UNSTAGGERED CASCADE 


ОШ ое ып intinite cascade of thin blades arranged as 
Сао an Pigure 5. Blade chord is 2b and blade spacing is 
Iwao Script zero referes to the reference blade and 
SuUbceript one refers to the adjacent blade in the positive 
ordinate (n) direction. The reference coordinate system is 
fixed at mid-chord with the abscissa (&) in the direction | 
of flow. | 

Two n coordinates j pi Oí od with ης originating оп 
the reference blade and n4 on the adjacent blade such that 
"S I η 


Now the following nondimensional variables are defined: 


E | 
X = ë (ІІІ-2) 
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Y - ,m-0,1 ΠΠ» 


(ІІІ-8) 


Ав the governing equation, the Teipel approximation is 
used, 


e 


фу Оф + Ky-0. (ІІ-55) 


УУ 


This equation must be solved subject to the following 
boundary equations: 

l. Zero disturbance upstream of the cascade. 

... πι τεπρεπςν οὐ the blades. 

Ρα η ῥεουιοιπ το. this latter requirement can be 
expressed mathematically as follows. Along the reference 


blade equation (II-50) may be written 


9 Ὦ | 
m z Ж O : | Ν 
9. (X0, Y570) та Ὁ. (III-5) 
and along the adjacent blade 
| oh, 
Y, OX, Y ΓΕ) =P O - IE O») 


To solve this boundary value problem, the potential 


function y is written as the summation 


O 
- O 
Еа Ос ку) + осу) 


+ 0, (X,Y7) (111-7) 
where ФО апа $T are the isolated blade solutions for the 
reference and adjacent blades and P and $. are interference 
potentials. 
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Ebenmonesopsoaech nas previously been proposed for the 
subsonic and supersonic cases by D. N. Gorelov [Refs. 3 and 
Bar Also refer to Samoylovich [Ref. 18]. 


Inserting (III-7) into (III-5) and (III-6) and realizing 


O 
ey 000) 7 v m-0,1 , 


m 3 


ойе obtains for the reference blade, 


= O ree e ME = 
'. UT 2, (Х,У КУЕ 21, 05 Y R) 0 
(III-8) 
ШШ бг the adjacent blade, i 
O = = = = 
ΛΈ R) * o ο DEL 2,5511 0) 0 
(111-9) 


The isolated blade solution is already known from 


Teipel [Ref. 9] as 2 ү2 
ос) 5:50 
ERY) pe τ (S)- Е ЕС ds 
| Е Š a a V —————— ———————— — 
E р σι ° /х-4 
Grate) 
f 2 
2 ү 
апа ΕχΡ[Σ. (χ- 5) - 2 τε] 
Ο al £ Q Е 
A (i) ————————— «s 
Υπο -ᾱλ VECES 
Са) 


The interference potentials must have the same form but 


contain unknown source distributions Uc and uj; ШЕШЕ, 


| 2 
2 
К Q O 
X pop c eee. x] 
o (X,Y )= - giis f u (s) MEE eor as 
9 va -1 ЕЕЕ 


(ІІІ-12) 
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2 
and 2 Y 
s > P Б EXP[“5(X-S) - E xu | 
NC --- "IDEM o X-S qs. (IIT-13) 
1 1 πο -1 21 τς 


Enese source distributions are approximated by the power 


series as follows: 


uU = D x Πα ιτ) 
апа 
Lo N 
Der (ІІІ-15) 
т N-0 mN 


Inserting these approximations, (III-14) and (III-15), 


into equations (III-8) and (III-9) results in 


L-1 N R fa x sN Ro 
) { DNA + 5 Ort Pin?) | ETETZ exe ESO -ϐ)- στα as 


0 Gobi 100) 


for the reference blade and 


1-1 [5 X N 2 2 
- 0 (III-17) 


termine adjacent blade. 


36 





—l .- 
— ο ας βαν» y | 


| πμ 
| — — ART Adjacent Blade 


LE mot 
Mena NEN wh you s 
т.т» Xo 
Reference Blade a 


Figure 6. BOUNDARY CONDITIONS FOR W PHASE LAG 


The coefficients in (III-15) may be evaluated for the 
case of simple ος, ЛО translation, ала mv phase lag 
between blades as follows. One arbitrarily sets the angle 
of maximum amplitude for the reference blade at an: = -1 
[Refer to Fig. 6]. Then the angle of the cent blade is 


a, = 1. For small a, equation (III-5) may be written 


V (X,Y = 0) = -1 - 1KX πι πο 


and (CD) 


0 (X,Y,= 0) = 1 + 3KX . (ІІІ-19) 


οἱ 





memee, the coefficients are 


- 


0 Ξ -] 
OO 
951 = -iK 
Ө (τος) 
Jo = 1 
and 041 = ik 


early, for in phase oscillation, the coefficients аге 


Me NO EE 


(ITI=21) 


and 91 = 0-1 = -jK . 


Satisfying boundary equations (III-16) and (III-17) at 
E points along the blades generates a set of 2L equations 
, νι. unknown interference coefficients. The coefficients 
ШІ Папа, 7 6,736) in (III-7) becomes a known function. 

computer program for the solution of this problem in 
the T phase shift case is developed in Appendix A. Some 
sample results are presented in Table II for the steady case 
and Table VII for the unsteady case. 
B. EVALUATION OF THE INTEGRALS OCCURRING 

IN EQUATIONS (III-16) and (III-17) 

The integrals occurring in (III-16) and (III-17) are of 


the form 


X oN Ke Q R2 | 
Iw = 4 (x-8)% EXP[—5(X-S) Е Y στὴ ds . (ІІІ-22) 


EyemMaking the substibution, О = XK-5, one obtains 


_ X*l(y ϱὸΝ 
N U 3⁄2 


2 2 
K QR 
J EXP[^g Ше n T dU | 009) 
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Recalling 9 = À + ƏiK, Iy is complex and must be resolved 


into real and imaginary parts for evaluation. 


The exponential term in (III-23) may be written 


2 2 2 2 
h Q R = еа _ AR 
EXP[^g U - т ч) ЕХРГ eo: Жү: 
. 3 2 
КІ KR 
“ЕХРІ E 5) (ІІІ-24) 


EJcEnow defines the following terms: 


Ke 
n 7 
À + HK 
"s 
y a tig 
2 


Equation (III-24) may be rewritten, 


2 2 
ЕХРГ 2-45 


4] = EXP[AAU - SEXP (i[-2KAU - ËJ} eee 


HO 


Separating (III-25) into its real and imaginary parts, one 


obtains the following two fundamental integrals from 


e T-23): 

Е u) 
Inn = S EXP[AAU - τ COS[2KAU 4 =] dU (III-26) 

0 
x+] 

Бито Quy EXPLAAU - B] SIN[2KAU * p] dU — (III-27) 
0 

where Jy = Ive = il (ІІІ-28) 


ВЕ 





Expanding (UN in the finite binomial series, 


ж- 


ШТІ 


М _ 
(Х-0) =) GT eye) (111-29) 
Ј= 1 
equations (III-26) and (III-27) become 


A un (13371 y; yN-J41 T NÉS ona. Ë. 
ΝΕ „2, EDT AED 5 Ü 


COS[2KAU + a au 
(IIT 30) 
and 


+ a N-J41 
ш Ч 


Tr s 
d N! X В 
бт é ΠΤ NIFI " fu с 


SIN[2KAU 4 a dU 


(III-31) 


The fundamental integrals in (III-30) and (III-31) do 
not lend themselves to analytical integration. Numerical 
integration techniques may be used but care must be exer- 
EM ans tnesveerrnity of the lower limit. Here the integrand 
is singular when J «5/2. The problem is further aggravated 
Metne rapidly oscillatory character of the functions 
SIN[AAU - Р] and COS[AAU - fJ. 

Three specific cases are posed for the treatment of 
these integrals. 

1. The Steady Case 

Setting K = 0 reduces (III-30) and (III-31) to the 


following: 
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IUE J-1 N-3+1 X*1 y 
NR ) "yr SW / να i [E 97 ау (ІІІ-32) 
J=1 А : 0 


ЕД 
l 


апа 


МІ 


Poe Oneri cal integration of this case is straight- 
τα. The exponential term forces the integrand to zero 
Maine origin permitting one to neglect that portion of 
the integral very near U = 0. Equation (III-32) may also be 
evaluated by expressing it in terms of the error function. 
(Note: complete agreement between the two integration 
techniques has been obtained thus verifying the numerical 
integration procedure for the steady case.) 

2, The Oscillatory Case with A = 0 
In the unsteady case for A = 0, equations (III-30) 


and (III-31) contain the basic integrals 


Х+1 Π.Δ E 
ΠΤ ση5οκλ + 7] dU (ШЕР 5) 
0 
X*l 9-2 sınfarau + EI au ΙΙ 
and U ` 
0 - 


Lora rotas ne previously mentioned difficulties near 
the lower limit, the integration must be accomplished іп 


WO Steps, 1.6., 
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ντι 


£ 
12-2 COS 


SIN 


SS 


E 
єт 2KAU 4 1140 + 


Р а 
F UC u 
0 0 


X+1 
Γι 
8 


J-2 COS 


E 


(1 2735) 


For judicious selection of (a), the second integral 
cae RAS may be numerically integrated in a straightforward 
manner. The integral on the limits 0 to (a), which will 
ιαπ τε referred to as the correction integral, may be 
approximated for sufficiently small (a) as follows. 


Recall the following relations: 


COS[2KAU + n = COS[2KAUJcCoS [J : SIN[2KAU SINE FI (III-36) 
and 
SIN[2KAU + =} = SINT2KAUJCOS[ E] + COS[2KAUJSINLZ] (111-37) 


and the truncated series expansions 


2 
COS[2KAU] =1 - n πι ο 
апа 
3 
SIN[2KAU] = 2KAU p. ατα ο 


Equations (III-38) and (III-39) are valid for 2KAU suffi- 


ciently close to zero. One may now write 
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Em É 
ie - COSE2KAU μου, 


0 
a 5 . 8 23 
| 972 соз| ав - (>кА)| U773 sIN[pJaU 
0 0 
2 а 4 
B (ҚА). |9972 сова 
0 
SE: 4 
6 0 
апа 
E 
f 0772 SIN[2KAU + ËJaU = 
0 
ἃ 125 а уз 
f 097 SIN[ËJ4U + (2KA) / 50772 совгЕ лас 
0 0 
(EN) NES T NE 
eo) USING lau 
0 
a E 
EC COS dU (ITI~41) 
πο 


Imzscheory integrals of the form 


a 
í ym/2 COS 
0 


Е 
stn 52140 


may be expressed in terms of the Fresnel integrals. Hence, 
integration of (III-40) and (III-41) can be accomplished in 
"closed" form. (In actuality there is a computer accuracy 


limit on this technique as described in Appendix A-IV.) 
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I S (l Oo of this technique for J = 1 was accomplished 


by checking it against a known analytical function as 








follows: 
a 
fee CCS EJ COS ji . 
j E gry KAU + Fld rea d Br + 180 = 
"UU a 
Ј —5- 12ка + 5Јаџ | (ENS OS 
ЕЕ 


where the correction integral was evaluated using the Fresnel 
integral technique, the second integral was evaluated using 
numerical integration techniques, and the integral on the 


Шы. Пав evaluated using the following relation: 


J mr «το[ῦ 1 π]40 - γπ/ο SOSIT * V2gb 1. πι ο 
For the case К - 0,01 and a = 0.01, Table ILA presents 
the comparison of the value of (III-43) with that obtained 
from (III-!2) for various R's. (The numerical integration 
was truncated at an Upper lame of 500 and the convergence 
criterion EPS for subroutine ZSUM was set at 1.0E-4.) 
Table IB displays the computed integral using (III-35) 
Пп ап Upper limit of X+1 = 2.0 for various values of (a). 
(The convergence limit placed on ZSUM was set at 1.0E-5.) 


Table variables are listed below. 


DK = K 
DL = À 
RHO = R 


hl 





Na Correction integral upper limit.) 


SCOR = Value of correction integral (III-41) for I = 1. 


CCOR = Value of correction integral (III-40) for J k: 


EN Value of numerical integration term from (a) to 
сее иеси оррет limit cor the СІМ integral. 


SUMR = Value of numerical integration term from (a) to 
see Upper πιπαις tor ne COS integral. 


SINT = Value of (III-35) SIN integral for specified upper 
Mimic and J =li 


CINT = Value of (III-35) COS integral for specified upper 
Eini crand J = и. 


SSIN = Value of (III-43) SIN integral. 


Eros = Value of (III-43) COS integral. 


3. The General Case, K Z 0 and À Z 0 


In the apparent absence of a suitable correction 
sc n q] τος this case, one must evaluate the full equations 
(III-26) and (III-27) by numerical methods. The difficulty 
close to the lower limit is alleviated to some extent by the 
ΙΙ опе а term. The exponential function drives the 
amplitude of the integrand to zero more rapidly than a 
negative power of U can amplify the highly oscillatory 
Memavior or the transcendental functions. For this reason, 
ο (ле case Of K = 0, one may justify neglecting that 
ποι the Integrals very near zero. 

nor resulosrare presented, the computer 
routine developed in Appendix A is designed to carry out 


Such a case. 





TABLE 


νεο τος ο CORRECTION INTEGRAL 


SCOR SUMI SINT SSN 
„5 34.062 1.016 33:03 Son OW 
- 14.45 3.458 ИЭТ 117285 
- 4.968 "τ 12.00 11.94 
2. е 10.16 9.045 8.987 
CCOR SUMR СІМТ CCOS 
0.5 15.48 19.80 Eun ms 25532 
1.9 πι 182 19.34. 17.56 17.60 
. ην 17.45 17265 11.69 
2. -4.490 КӘНІ? 8. 696 8. 736 
ο ο ο ος ναι ος» COMPUTED USING 
USING (111-42) FOR AN UPPER LIMIT OF 500. 
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И ЕКЕШ ИЕБИ ЕК ТЕТСАТТОМ OF CORRECTION INTEGRAL 


DK= 0.500 DE 9-0 RHO= 0.500 Al= 0.010 
SCOR= 1.57005D 00 CCOR= 1.671140-01 
SUMI- 3.968510 00 SUMR= 3.240130 00 
SINT= 5.538560 00 CINT= 3.407240 00 
DK= 0.500 μμ. 5 ΒΗΟΞ 0.500 AL= 0.050 
SCOR= 1.66484D 00 CCOR- -2.65483D 00 
SUMI: 3.873720 00 SUMR= 6.062070 00 
SINT- 5.538560 00 CINT= 3.407220 00 
DK= 0.500 Ш- 000 RHO- 0.5009 ο ο 
SCOR= 3.74310D 00 CCOR= -1.12801D 00 
SUMI- 1.79545D 00 SUMR= 4.535260 00 
SINT- 5, 538560 00 СІМТ- 3.407240 00 
DK= 0.500 DE 2020 RHO= 0.500 AL- 0.200 
SCOR= 4.620260 00 CCOR- 4.962700-01 
SUMI= 9.182970-01 SUMR= 2.910970 00 
SINT- 5.538560 00 CINT= 3.40724D 00 
DK= 0.500 πο RHO= 0.500 ME 002300 
SCOR= 4.87826D 00 CCOR= 1.27508D 00 
SUMI= 6.602940-01 SUMR= 2.132160 00 
SINT- 5.538560 00 OT 245500 
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See CALCULATION OF THE POTENTIAL FUNCTION AND 
THE PRESSURE COEFFICIENT 


The potential function is obtained from equation (III-7), 


MES “9 
V(X,Y,) = De + De, 4 $4 (X,Y1) 1 $4 CX,Y4) 
(III-T) 


The interference potentials Ὁ and %. пау now be 


1 
computed using the previously determined coefficients (DN? 
DEEENMUU N- 0,..L-l) in equations (III-12) and (III-13). 
Beustjion (III-T) may be rewritten in terms of the coeffi- 


cients as follows: 





2 
| L-1 X “N 2 Y 
v(X,Y.)= J \- BEE pe Ges). 0 Jas 
ο "τ. е Q 1 (Х-5) 
2 
X <N 2 
1 5 K Q 1 
+ E (95 yt? a? J I ЛЕРДЕ 5(Х-5)-Т нуге! 
(ІІІ-!!) 


Тһе resulting integrals may be evaluated using the same 
aenn dues utilized in the coefficient derivation. 


The pressure coefficient is obtained from the. relation 
Cp - x LEUTE (III-45) 


Ву placing the partial with respect to X of equation (III-44) 
into (III-45), one may obtain the pressure distribution at 
any point in the flow field as well as on the blade surfaces. 
Further, the force and moment coefficients are easily 


obtained by integration over the blade chord. 
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E VACUATTON OF THE STEADY STATE 
POTENTIAL ALONG THE- BLADE 


MancEcenpucation ofi the steady state potential along the 
Uer surface of the reference blade follows from equation 
(I11-44) evaluated at K = O, БЕ ТТ Ci nie 


development which follows is for the m arse seine ea ae 


where [From Eqns. (III-20).] 





aoe al 
OO 
ή, ΤΙ, 
θη ο = d 
Equation (III-44) becomes 
2 
\ Β 
_ 1 πη Х ЕХР{- р =] 
Vena (X50) = — —— 05 + f eee е 
ЖІ zn voce 
ee αυ 2 
1 5 R 
+ — P D f —— EXP[- 7$ v—JàS 
ТТ κο IN 1 y = 
| oN | 
- ——— dS (III-47) 
жу с 


The computer program is documented in Appendix B-I. The 
resulting potential values are compared with those obtained 
Ет спе Sandeman wall interference solution of Section III-E 


шш Тар1еес TIT and IV. 
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E. LAPLACE TRANSFORM SOLUTION FOR THE SONIC WIND 
TUNNEL WALL INTERFERENCE PROBLEM 


A parallel development to the collocation technique is 
mona in the Laplace transform solution of the sonic wind 


tunnel wall interference problem [Ref. 4]. 


= sam 


= Ба a — po == Su Virtual Blade 


R = 2d 





Tunnel 
Wes : ο Walls 





Figure 7. SONIC WIND TUNNEL 


Figure Т shows the arrangement of a blade facing a solid 
wall. Due to symmetry, the solution for the flow between 
the blade and the wall will be the same as that which would 
exist between the blade and the virtual blade beyond the 
wall. Hence it is the same solution expected of an unstag- 
gered cascade with π phase shift between the blades. 

The basic equation for oscillatory sonic flow from 


section II is 2 
Pes - Оу, + Ku = 0. το) 
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The boundary-value problem can be formulated as follows: 
ШӘ ТИІСЕ is no disturbance upstream of the blade 
leading edge. 
2. Flow tangency race must be satisfied along 
pe blade and along the solid wall. 


too latter condition may be stated mathematically, 


-v (X) at Y 


N 
© 


(ІІІ-48) 


Y, (X,0) 


1 
= 


and Y, (X,d) at Y = sd (БІНЕ ЛО) 


Denoting the Laplace transformation of a function with 


11ο X by an overbar, i.e., 


со 


f(p) = L TX) = f εν πο αχ. 
0 


one obtains for the Laplace transform of Eqn. (II-49) 


== а ! 
- À“) = 0 ІІІ-50) 
Vay Y ( Э 
Ec 2 
where М" = др - K (ІІІ-51) 
mie boundary conditions transform into 
I s Y = 0 (ІІІ-52) 
апа b; = 0 Y = +4 (ІІІ-53) 


ΠΠ το το) засів пуз пе (111-52) апа (111-53) 


has been given by Miles [Ref. 19], 


Ф(р,Ү) = = ¥(p) CSCH(Ad)COSH[A(d-|¥| ) Isen(¥) НЕСІ 
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ΠΕ oneself to the flow between the upper tunnel 
ЖС = A and the blade one has Y > 0. Hence the absolute 
Value sign may be dropped. 

Mommnavert (111-54), a suitable expansion of the equation 
Was given by Miles as 


со 5. 


EO Y) = Y v(p) ) (EXP[-A(2Nà4Y)] + EXP{-A[(2N+2)da-Y]}) 
N=0 | 
(III-55) 





mecalling the Laplace transform of 
eS y 1? of Ix Rea > 0 
being 
= 2 уар 
Т (р) = | е P^etx) dX = Ge Re pp > 0, (ІІІ-56) 
0 | үр 
ame has for 
F(p,Y) = 5 EXPL-A(2Na+Y)] 
e 2 
where À = yQp-K 
the inversion κ Q (2Nd+Y)® 
| ЕХР(-< Х p] 
Y) (ІІІ-57) 
TQX 


mietretore, with the aid of the convolution theorem, one 


obtains the following solution: 





ШУ) = T г ey ExPLÉ (х-8)] 
| үт © " y X-S 9 
со 2! | 
ο Cer + exer. g Legen) as 


(III-58) 


e 





Пс blade surface where Y = 0, this result simplifies to 








| qa 
fee) = - — ue ЕХР[ 7; (Х-8)] 
Итә -1 X-S 
= оқа τ 
) 4 EXPL- S—a—] * EXP[- xc i 
а Х-5 T 
N= 0 
(ІІІ-59) 
For steady flow where K = 0, one obtains the further 
Splat ication, 
Х 
ПЕР) - - Е | V(3) as 
Үт ~ ] Х-5 ы 
со x p 
"τι; NS ав }. (III-60) 
N=] zo ve 


This expression was previously given by Sandeman [Ref. 4] 
in his analysis of choked wind tunnel flow. Equation 
(III-58) represents the generalization to oscillatory flow. 

Appendix B-II describes the translation of (III-60) into 
computer program. The solution of this function provides 
press teady state check of the collocation technique as is 


Shown in Tables III and IV. 
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IV. PRESENTATION OF RESULTS 


fable Ll presents the interference coefficients for 
the steady state case obtained for the controlling param- 
eters indicated in the tables with the computer routine 
developed in Appendix A. 

Table III contains the value of the perturbation 
ШЕСЗепгігзі along the reference blade obtained from the 
routine developed in Appendix B-I and utilizing the inter- 
ference coefficients of Table II. These results may be 
compared with those in Table IV obtained from the steady 
meee Laplace transform solution for the sonic wind tunnel 
wall interference problem of Section III-E and programmed 
in Appendix B-II. 

Table V demonstrates the erroneous correction integrals 
obtained from Subroutine COR of Appendix A as the result 
Eu ccmpoter numerical inaccuracies. A comparison of the 
interference coefficients for two cases before and after 
the anomaly was detected is presented in Table VI. 

A sample of the unsteady interference coefficients "ος 
the case À = O is presented in Table VII. 

Table VIII compares a sample steady state output with 
the output of the general case À Z O, K Z 0, for very 


Small K. 
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Bement TAC SOLUTION FOR THE FOCCOWING PARAMETERS 


DL = 0200 
THE 

ШЕР (00) --2 
ШЕГІГОІ)--1 
ШКЕ {О2)= 3 
ШЕР (905) --1 
DEL(04)= 1 
MEE TOS) =-6 
DEL (06) =- 2 
DEL(O7)= 3 
DEL(08)= 1 
DEL(09)=-9 


° е ө ө ^ e 6 ғ . 6 
-UU0)0- Ou Ono 
С оллоо Ce OU 
OO VAOU ORON 
JI МОСОСООСОФОУ, 
MODH Ae Mmm mm 
O CO (N CA une -e i 


DK= 


MEME 
ОЭООСООООО 


WMDUWONDUIAYO 
АЛЕ ΜΗ ων 
| 
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' ΠΓΕΕΒΠΕΕ COEFFICIENTS 


мш менн оо 


RHO= 1.000 
ΕΙ (10) = 2.0453110E 00 
DEL(11)= 1.2504177E 00 
DEL (12)=-3.0618846£-01 
DEL(13)- 1.5100145E-01 
DEL(14)=-1.0638958E-01 
DEL(15)= 6.7890167E-02 
DEL(16)= 2.0572830E-02 
DEL (17) =-3.3797889E-03 
DEL(18)=-1.3067639E-01 
DEL(19)= 9.7604275E-02 


ШЕШЕ ОБ THE POTENTIAL ALONG THE BLADE 


STAT ION 
Х--0.90 
А-- 0.60 
Х=-0. 70 
Х--0.60 
Х-- 0.50 
Х--0.40 
Х--0.30 
Х--0.20 
Х--0.10 
X=—0. 00 
X= 0.10 
Х- 0.20 
Х- 0.30 

- 0.40 

= 0.20 
X- 0.60 
Х- 0.70 
Х- 0.80 
Х- 0.90 
Х- 1.00 


p 
POT= 
POT= 
РОТ- 
РОТ- 
РОТ- 
POT= 
POT= 
POT= 
РОТ- 
РОТ- 
РОТ- 
РОТ- 
РОТ- 
РОТ- 
РОТ- 
POT= 
POT= 
POT= 
POT= 
POT= 


ШЕТІ 
12173099 Е 
DAMM MOMO 32 E 
ESA 
4.1155638Е 
91750020F 
O ZE 
ποσα Ε 
8.1755124Е 
ΠΕ 
ІШ ο) ο {ΓΕ 
ISA 
1.211 TOE 
Io» TOf 
Pettis DOLE 
1.51 49505 
1.6175568E 
J./1715552E 
1.81754#T6E 
19153308F 
РЛЫ ФЕ 


TABLE III 
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00 
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00 
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01 
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01 
01 
01 
01 
01 
01 
01 
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ο να SOLUTION FOR THE FOLLOWING PARAMETERS 


DL- 0.200 


DK 


NJ NI UERSFERENCE COEFFICIENTS 


CO O ONA. IO. O 
ло Алом ОМ 
el 3 OA! OO RO U1 PO 
WRIWNUIO OO ш) чә 
INWWINMAJaIN — O — 
Mm UTU UIN) — J vO O 00 
Ολο о. ы 


© о e 9* e ° o е o 


EE COO) -—-8 
DEL(01)=-6 
DEL (02)= 1 
DERCOS) =-6 
DEL(O4)= 2 
Milla) 5 ) =— 1 
ШЕ! (26)- 5 
МЕЗО /)=-5 
DEL(08)=~4 
DEL(09)= 5 


Mmm?mmmmmmm 
! 

ο ον S G Ὁ στο © 

NOMI NO RO RO AOI нә не ке 


= 0.0 RHO= 2.000 
DEL(10)= 8.0229187E-01 
DEL(11)= 6.1007094E-01 
DEL (12) =-1.4253497E-01 
DEL (13)= 6.6098273£-02 
DEL(14)=-2.7520720E-02 
DEL(15)= 1.2785759£-02 
DEL (16)=-5.6872237E-02 
DEL(17)= 5. 7169359E-02 
DEL(18)= 4.0132940E-02 
DEL (19)=-5.0579518E-02 


peeve OF THE POTENTIAL ALONG THE BLADE 


STATION 
Х--0.90 
0.50 
х=— 0.170 
Х--0.60 
Х--0.50 
Х--0.40 
A0, 30 
Х--0.20 
Х--0.10 
Х--0.00 
Х- 0.10 
X= 0.20 
X= 0.30 
Х- 0.40 
Х- 0.50 
Х- 0.60 
Х- 0.10 
Х- 0.80 
Х- 0.90 
Х- 1.00 


POTENTIAL 
РОТ- 8.1581414E-01 
POT= 1.3204842E 00 
POT= 1.8202248E 00 
POT= 2.3200750E 00 
POT= 2.8200855E 00 
РОТ- 3.3201075£ 00 
POT- 3.8201075E 00 
POT= 4.3200874E 00 
POT= 4.8200750E 00 
POT= 5.3200960E 00 
РОТ- 5.8200798E 00 
POT= 6.3200665E 00 
POT= 6.8200712E 00 
POT= 7.3200731E 00 
POT- 7.8200769E OO _ 
POT- 8.3200636E 00 
РОТ- 8.8200598Е 00 
РОТ- 9.3201199Е 00 
POT= 9.8202019E 00 
POT= 1.0319364E 01 
Ῥ ο ΤΙ 
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POTENTIAL S OLUT ION FOR THE FOLLOWING PARAMETERS 


DL= 0.200 


Ше ГГТЕКРЕКЕКМСЕ COEFFICIENTS 


ΕΦ СММ мнн 


l l H H H H dg dg gu H 


Са са a a e 
mmmmnnrmmmnm 
[s = r= (ua ЕКТІ. 
қ ут та жа /-- s s w j 8 р 
СОС ο S 
Осо !С Un Un To 


NIT) Co UND m 00 ps 
Μου O awl WN 
FUI O US UI Ly 
SOI Ore fb UI 
N.N Cs Q) çp + O + 
лм рх Уо Ыл 


€ 0 9 9 е ее еее 


DK- 0.0 RHO= 5.000 

SE-01 DEL(10)= 1.1538458£-01 
6 E-01 DEL (11)= 1.8925136E-01 
6E-03 DEL(12)-—7.1354136E£-03 
5E-02 DEL(13) 2-2. 2944096E- 02 
9E-02 DEL(14)= 2.5751479E-02 
15-02 DEL(15) =-2. 870987 1E-02 
1E-02 DEL(16)= 6.209834 7E-02 
1Е-92 DEL(17)=-4.1337617£-02 
OE- 02 DEL(18) =-5. 2630220E-02 
9Е-02 DEL(19)= 5.1547259E-02 


MAGUE OF THE POTENTIAL ALONG THE BLADE 


STATION 
X=— 0. 90 
X=—0.60 
“τος {ο 
ΧΞ- 0.60 
Х--0.50 
X=-0. 40 
X=—0. 30 
Х--0.20 
Х--0.10 
Х--0.09 
Х- 0.10 
Х- 0.20 
Х- 0.30. 
Х- 0.40 
Х- 0.50 
Х- 0.60 
Х- 0.70 
X= 0.80 
X= 0.90 
X= 1.00 


POTENTIAL 


DONE 726253€=0] 
РОТ- 1.1282368Е 00 
POT= 1.3855448E 00 
РП = ӘСЕ ОЕ οῦ 
ος ο ορ ο E ОО 
POT ο οσοι ας ОО 
РОТ- 2.2310429Е 00 
РОТ- 2.4321480Е 00 
РОТ- 246326466 00 
POT= 2.8328724E O0 
РОТ- 3.0329742Е 00 
РОТ- 3.2330236Е 00 
РОТ- 3.4330444Е 00 
POT= 2.6330500E 00 
РОТ- 3.8330574Е ОО 
РОТ- 4.0330648Е 00 
РОТ- 4.2330627E 00 
РОТ- 4.4330502Е 00 
РОТ- 4.6330080Е 00 
POT= 4.8337231E 00 
TABCE IIT 
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mare tat SOLUTION 
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Perth FOLLOWING. PARAMETERS 


RHO= 10.000 
DEL(10)= 1.5652659E-03 
ПЕ (11) = 8. 5002407203 
DEL(12)= 1.4888730E-02 
DEL(13)= 5.3168200£-03 
ΡΕΙ (14) -- 7.090728 7E - 03 
DEL(15)= 6.0283695E-04 
DEL(16)= 4. 9886145£-03 
DEL (17)=-3.8378909E-03 
DEL (18) =-1.6985249E-03 
DEL(19) = 2.1847319E-03 


Mabe OF TRE POTENTIAL ALONG THE BLADE 


STATION 
X=—0.90 
Х-- 0.80 
ο. ο 
ΧΞ-0. 60 
5-50 
Х--0.40 
Х--0.30 
Х--0.20 
Х--0.10 
Х--0.00 
X- 0.10 
Х- 0.20 
Х- 0.30 
Х- 0.40 
Х- 0.50 
X= 0.60 
X= 0.70 
X= 606.80 
X= 0.90 
X= 1.00 


POTENTIAL 

POT= 7.9786545E 
POT- 1.1283703E 
POT= 1.3819685E 
POT= 1.5957613E 
POT= 1.7841234E 
POT- 1.9544506E 
DO E 
POT= 2.2573242E 
РОТ- 2.39499096 
POT= 2.5258138E 
POT= 2.6510839E 
POT- 2.7718019E 
POT- 2.8888006E 
POT- 2.0027514E 
POT- 3.1142006E 
POT= 3.2235994E 
POT= 3.3313122€ 
POT= 3.4376440E 
POT= 3.5428400E 
POT= 3.6471376E 

TABLE III 
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ο δα =10d 0670 =X 
ος πε ο = 10d OT°O =X 
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V. DISCUSSION OF RESULTS 


A computer program was designed to compute interference 
coefficients for the collocation technique developed in 
Section III-A. For the steady case, the resulting coeffi- 
cients in Table II were applied to the technique developed 
section III-C and III-D to provide the values of the 
Busco surface potentials in Table III. For this same case, 
tne sonic wind tunnel wall interference solution of Section 
III-E provided the data presented in Table IV which verified 
МЕ Исо осатсіоп technique for that case. 

Although no "separate source" checks were made for the 
coefficients presented in Table VII (The case of simple 
oscillation about the mid~chord with no translation.), every 
effort was made to verify the procedure. Table VIII com- 
pares the proven steady state solution with a general small 
σσ Unsteady integral approximation techniques were 
Jüstified in Sections III-B-1, -2, and -3 as summarized 
De ION. 

Premreulbies arising from the solution of the occurring 
integrals in the unsteady solution were not insurmountable 


- 


but required careful handling. Whenever numerical integra- 


| : 5 J: SINE _ : 
tion of functions involving s ΟΥ Gos (5) near X = 0 is 
ΜΠ. cConpacy is Suspect. One must closely inspect 
ο ποιο limits, and results to ensure validity. For 


example, for the unsteady case and A # 0, with no apparent 
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means of computing the integral in the immediate vicinity 
Пе origin, one πη. τας, r agains. cHeWsmertcomines 
En computer capacity and programming methods, just what 
approximation must be applied and, once applied, how 
 гассе it is. Іп this case, the question arose as to 
ШЕ С105еїу one might approach the origin with numerical 
ШРерпгасіоп techniques before accuracy was degraded. At 
m απο Cime the required computer time increased 
dramatically as the integrand became more ill-behaved. 
Most difficult question that arose in this investi- 
ου "μας the accuracy of the correction integral in the 
Бе 1 = 0. Once this approach was found valid for J = 1 
ШЕ ОШ the calculations shown in Seetion III-B-2, there 
КОК ne need to improve the accuracy of the technique 
ENNCUSUure its validity for J > 1. In question were the 
Bemeer of terms needed in the truncated series shown in 
Equations (I11-40) and (111-41), the validity of the inte- 
Ma eduction formula used in subroutine RED, and the 
EENection of the upper limit of the correction integral. 
Duhouph four terms were used in the truncated series 
Шо Ое case J = 1 (Refer to Appendix A-IV), only three 
memmowaitected the significant figures of the results for 
Pie parameters used in this investigation. For J = 2, 3, 
ena 4, the number of terms used was decreased by- one for 
each computation. The primary motivation for this step was 
the inaccuracy of the integral reduction formula for large 


шанс агре К. 
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Referring to Tables VA and VB, it was determined that 
uu col inaccuracies began to occur for large J in the 
En cEBRSL reduction routine. The effect of this inaccuracy 
was not important for small R but for the case of R = 5, 
for example, it overwhelmed the results. For integrals 
Se che form 

AL 


/ U 
0 


T ST] 
J-£ SIN 


cog (E/U) au, 


becomes large and for AL << 1, the value of the integral 
Kurs come small. However, in the case of R = 0.5, the 
magnítude of the integral values diminished with increasing 
"uto J = 9 but increased for 7 > 9. Іп the case of 
OS. the anomaly in the computation affected integrals 
l ас low as five. (In this latter case, it is obvious 
ші пе erroneous answer for J = 13 must affect the results.) 
The cause of this error was determined to be differences 
taken between large numbers in the reduction formula. The 
ΜΙ obtained for the Fresnel integrals was only as accurate 
И согуе fit used in Subroutine CS. This error and the 
Pe@omeeeneraced in computing the SIN and COS functions 
υπ πο πρ in the formula resulted in erroneous answers of 
wn eaa cude that was directly related to J and R. To 
τοσο the problem, the use of the integral reduction 
routine for J > 5 was eliminated. 

ШИЕ, пппегіса! techniques were resorted to. Неге, 


MEDIE the amplitude of the integrand was small 
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permitting one to neglect that portion of the integral near 
zero. This technique seemed to produce reasonable results, 
at least the results did not have large magnitudes that 
ООЗЕ have affected the final answers. The resulting 
integral values were in some cases less than the convergence 
limits placed on the numerical integration subroutine, ZSUM, 
EBnoEwere, therefore, assumed to contain error not greater 
pen che said limits. Further, as these integral values 
became small, it should have been possible to Ne eee chem 
EXMUosetchepr since the values of those integrals with small 
(J) predominated in the series. 

Comparative results for the case K = 0.5 and R = 0.5 
and 5.0 are presented before and after the error was 
ου in Tables VIA and VIB. 

The upper limit of the correction integral was selected 
such that the approximation 92s (E/U) = Sos (DU + E/U) is 
ИЕСИ the range 0 < U < AL. The point AL also needed to 
πια the capabilities of the numerical integration 
ο πε οῃ the limits AL < U< X + 1. Such a point was 
picked at AL - H/100 which ensures that (For U « AL) E/U 
was at least 1.0E4 larger than DU. Such an AL was also 
nin ohne Capabilities of the numerical technique over the 
. προ Of parameters used in this investigation. Table IB 
ze monstraves that the use of the truncated series in this 
approximation decreases its sensitivity to AL. 

Computer time could have been reduced with a larger AL 


if one were willing to accept reasonable inaccuracies as 
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ШОС еги DU became more significant. A smaller AL increased 
Fomputer time апа the πα πας ОЕЕО numerical integration 
technique became more and more questionable as the integral 
meeame more ill-behaved. Computer time for a single run 
the program (for à = 0) on the IBM 360/67 computer 
stem was on the order of 15 minutes. 

It should be noted that in the absence of a good error 
 /гпгепепс yardstick, statements concerning more or less 
Aita re indeed relative. It is difficult to assess 
Dur-ccomputer numerical error or approximation error be- 
comes more important. It is for this reason that an inde- 
E u-chneck of these results against another source is 
Necessary. The Laplace transform solution for the sonic 
unsteady wind tunnel wall interference problem would, in 
ου s y be such a check; however, this technique will 
result in integrals of the same form as those used in the 
collocation solution and will require similar (or the same) 
on techniques. It will be difficult, then, to obtain 
КИ елу independent check of these results. 

Once verification οι πο κο onto results presented 
ου τς investigation is in hand, application of Spreiter's 
local linearization technique alluded to in Section II 
should be applicable and lead to reasonable analytic 
results for transonic pressure and moment distribution 


prediction. 
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ҮІ. SUMMARY 


Following the basic approach offered by Gorelov for 
subsonic and low supersonic flows, a collocation technique 
was applied to the purely transonic case superimposing the 
isolated foil solution provided by Teipel with unknown 
interference potentials. A computer: routine was developed 
which woulá generate the coefficients for the unknown 
Iu -rrerence source distributions in the steady case. The 
 ГегГегепсе coefficients in hand, the potential along the 
Dlade was determined through a second routine. These 
ΙΙ; were favorably compared with the potential solution 
derived from a routine based on the Laplace transform 
treatment of the sonic wind tunnel wall interference problem 
for the steady case. All three of these computer routines 
initially utilized series expansions which reduced the 
Seeurring integrals to functions of the well known error 
ΜΝ ου COmlirmation of these results led to an extension 
S e collocation approach into the unsteady regime. 

First attempts to utilize series approximations for 
Do eunt eady case resulted in failure due to nonconvergence 
of the series. Faced with a lack of an apparently simple 
ОШ "ο. τους resulting integrals, numerical integration 
was resorted to. Неге difficulties arose centering around 
Mer p=aunologlcal behavior of the integrands near the origin 


ОООО Л Ое ОП in the Discussion of Results. 
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As a basic numerical integration technique, the "pyramided" 
ENupson s rule method vas incorporated. For the steady case 
where the integrals were easily computed, the results of 
this method were checked against the steady state case 
already in hand and found to agree perfectly. 

li corestest difficulty occurred in the unsteady case 
for A = O0 where one could not reasonably neglect the inte- 
grals in the interval near zero. To overcome thís problem 
an analytical approximation (for the integral evaluated 
near zero) was introduced and checked against a known 
analytical solution shovn in Section III. The technique 
uu ved the reduction of the integrals to functions of 
EM nel integrals. However, this procedure resulted in 
additional probiems which were not immediately detected. 
ИШЛЕ the technique was theoretically sound, the computer 
Memetrical error pointed out in the Discussion of Results 
Meecame apparent after study of the output. Once the cause 
was detected, the anomaly was remedied in a simple manner. 

The της for the unsteady interference functions 
as presented here shall require a "separate source" check. 
Ше ТШартасе transform solution for the oscillatory transonic 
wind tunnel interference problem should provide such a 
check. If these techniques verify one another, a third 
independent check could be the Hamamoto linearized high 
frequency solution. 

Development of required computer routines for the 


potential and pressure coefficients from the collocation 


85 


ы. X 





 ШБартасе techniques should be facilitated by the use of 
ENDPoutines derived for tnemünstreagyzesefficTjent routine 
presented heroe. 

Шо τοις that the collocation technique will prove а 
111 tool in transonic flow studies in that the free- 
dom which one has in specifying boundary conditions will 
ΙΙΙ treatment of problems which defy solution by other 
ИКЕ ез. For example, the sonic wind tunnel wall inter- 
ference technique is only good in the t phase shift case 
агат соттссасіоп technique will provide results for 
rsi se relationship. Further, it should be valid for 
Eo. ceomplex boundary eonditions sueh as those posed by 


pcunnels with porous walls. 
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APPENDIX А 


FORMULATION OF THE COLLOCATION COEFFICIENT ROUTINE 
A-I. MAIN PROGRAM 


ies purpose of the main program is the solution of 
Equations (III-16) and (III-17) for the unknown inter- 


ference function coefficients. 


N 2 2 


L-1 X 
N RIG S A: p 
KR. M A 7 s eger р 


(ІІІ-16) 


Bus 8 K^ πμ... 
ᾱ [δηιρ ΟΝ ΙΙ Е Јаз} =0 


(111-17) 


πο τος (III-16) must be satisfied for each station along 
the lower blade while simultaneously equation (III-17) must 
be satisfied at each station along the upper blade. For 
puecpurposes of this ιτ, е лос ауе been selected 
(X, a O 22. r = O,l,2,...,9). Computations near the 

Ἱ ο πα edge are avoided due to the singularity there. 


Пеи січеп К, A, and R, the term 


Ker eh 2 2 
RIS S K Q R 
ἘΞ f —À—MuEXP[-—(X-S)-45 w—]àdS = SR(X,N)+iSM(X,N) 
217 1 (xs)? Q 1 2-5 


ШЕГІН 
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Mea complex function of the parameters X and N. Тһе 
pmp lex value of this function at each (X,N) is provided 
by Subroutine WSUM and will be dealt with under Appendix A-II. 
WSUM returns the value of (A-I-1) as real and imaginary 
Meets in arrays SR and SM respectively. 

пе problem then is to evaluate the various known or 
computed values and insert them in БОЕВОЕ sultanes 
ο Όττο Simultaneous solution. 


erzerve (III-16) and (III-17) become 


9 
Σ (ы * (o A CSROGN) + ism(x,N))} = 
Mee QC 


and 


[DA КОО ОО О 4 ISM(X,N))] = 0 
(A-I-3) 


1 M\O 


N=0 


With 9 and 6 determined Dy boundary conditions at each 


ІМ ОМ 
of the ten stations, (A-I-2) and (A-I-3) form a series of 
twenty equations in unknowns Dow ΠΠ Ν᾽ Ша O 
fete, Lor the unsteady case, Ом and Zw aro complex. 

Taking the set of coefficients, Oyy» from терр 


~ 


(A-I-2) and (A-I-3) are rewritten for the in phase case as 


N . _ 2 
Pox + Day ( SR(X,N)+4SM(X5N)) | πι) DE МОС) 


I mn 


N 


+ i(SM(X,0) + DK SR(X,1)) τις 
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апа 


9 E 
i (αν δα ica] = (SR(X,0)-DK SM(X,1)) 


N=0 


Zst) + Dkwok( xX, 1) Caio 


where DK replaces K in the computer routine. 
For the vr phase shift case, the coefficients (III-20) 
apply and the lower blade boundary equation (A-I-2) 


becomes 


UDIO 


(D pX +2, (SR(X,N)+1sM(X,N))] = -(SR(X,0)-DK SM(X,1)) 
N=0 | 


- i(SM(X,0) * DK SR(X,1)). (A-1-6) 


The system of twenty equations may be written in matrix 


form 


μα] 

ol 
It 
ση 


CRUS 


where X, D, and B are the matrices of coefficients, unknown 
D's, and right hand sides. X, D, and B are complex and may 


be written 


οἱ 
| 
σσ] 
-- 
>l 
H 


η 2 
P = D. T Doi | ο 00) 
апа B = 81 + Boi. 





Men in (A-1-7) 


- 


or (D. X.-D x Wen ЫП Игр Е (A-1-9) 


or 


K =y D B 
(5 " | 2) B (>) Т 
X UBND B 


Jhe twenty complex equations in twenty unknowns 
represented by (A-1-7) have been transformed into 10 real 
equations in 40 real unknowns. One has only to solve 


mor D4 and D. and set 


D = D4 + 1D5. 


PROGRAM STEP ANALYSIS 
(Refer to Routine MAIN at the end of this Appendix.) 


Double precision is used throughout the routine for 
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merical accuracy. else no ο ο vario ls arrays and 
common storage areas that are declared will be detailed 
below. 

Various arrays ana fixed parameters are initialized 
in MAIN to minimize computational redundancy in later pro- 
πα. The arrays COC(10), COS (10), ΠΠ, 
eed in DO LOOP 100, are utilized in Subroutine WSUM. 
S DK = K, and RHO = R, the controlling parameters, 


are initialized and PI = T is set. The value for 


0.0 


A = Ne U TIO E 
(DL)? + 4(DK)? 


(А-І-11) 
Nx onputed. For the case DK = 0, A is set to zero to 
οποια possible computer problems involving division by 


ο ου Additional parameters are computed as follows: 
£ 


2 
T DL(REO) 

В = 2 ντ πο 
| DL(DK)* 

тШс лс у, SS τη, 15) 


Γι.» ; πι. 


5 
a 2 (DK) _ (A-I-14) 


(DL) + KR)“ 


2 
= = p (A-I-15) 
_ 1 2DK рес 
G = 5 TAN (FT (А-І-16) 
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Here G is set to 


G = 7/4 (S 


one case DL = 0 to avoid division by zero. 


AK = SIN(G) Cee 
BK = COS(G) | (1-1-18) 
2 „27% 
; . RHO LL: ! (DK)^] (Are) 
2 VPI 


For future ease in changing the number of stations on 
ШОШ аде, the parameter L is used. Here L=10is set. 


Ше іпісіг! station, X. = -0.9, and the interval between 


1 
Pwavtons, XINT = 0.2, are set. 

Miesiirst write statement is included to indicate the 
EsmUprolling parameters DL, DK, and RHO. 

Finally, the matrix manipulation section begins. 
Referencing (A-I-10), array AA(40,40) is dimensioned to 
hold the coefficient matrix. BB(40) holds the matrix of 
right hand sides before solution and the matrix of unknowns 
after solution.  (A-I-10) reflects the diagonal duplication 
of the real and imaginary submatrices, X and X, . Note 
that equations (A-I-4) and (A-I-5) form adjacent rows in 


the array and that there is duplication between them for a 


given X. The first two rows of X are as follows: 
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1 X n os SR(X,,0) SRQGQ,1) ... SR(X,, 9) 


l 1. 
© 
SR(X, ,0) SR(X, 51) T SR(X, 59) 1 X en X, 
(A-I-20) 
Boerresponding to the unknowns 
Leo дро1 E Prog Pao Dan u Das 


where Dun = Davy f ig: The first two rows of X, are: 


0 0 n 0 ο ο το С) 
ED) SICK, s2) ... SI(X,,9) 0 0 ee 0 
ОЛОТ; 


vorresponding to 


P: oo τοι 5 Pros Pin Pu: ps Pg 


Reference to (A-I-10) shows the placement of these rows. 


The first row of X, "ΙΡ ο ο ο. S well as 


EE n j iro AA(21,40). The first row of X. fills AA(1,21) 


2 
ο πο, and, with signs reversed, AA(21,1) to AA(21,20). 
For the T phase shift case, the right hand side matrix is 


filled four rows at a time as follows: 
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- 


ΕΡΤ) 


-(SR(X,,0) - DK SI(X,,1)) 


BB(2) = (SR(X, ,0) “Dr SI(X, ,1)) 
BB(21) = -(GI(X, ,0) + DK SR(X, ,1)) 
BB(22) = (SI(X, ,0) + DK SR(X, ,1)) 
and (А-І-22) 
BB(3) = - (SR(X,,0) - DK SI(X,,1)) 
p = (SR(X,,0) - DK SI(X,,1)) 
BB(23) = -(SI(X,,0) + DK SR(X,,1)) 
BB(24) = (SI(X,,0) + DK SR(X,,1)) 


ӘШШ 50 forth. 
EOZEEOOP 1000 on J = 1l,b-varies the station from X, to 


X It sets two rows in submatrices X, and Χο, ог, 


wor 1 

equavalentiy, four rows in AA and four in BB. A is used 
to set the power of X in submatrix Xi απ ο Ισ οὐ cet 
to A = xo = 1, JJ is a counter used to «ον ος ενος 295v 
Hu Succ 200- 1,3....,109-ав5 2 - 1,2,...,10. 

With Xy set at X = -0.9, WSUM is απ. το provide the 
integral term in the boundary equations. WSUM computes 
the integral for X4 πη. το: Ου ,..;: These ten 
values are returned in common arrays SR and SM. SR holds 
the real parts and SM the imaginary as shown in functional 


Bos ο SH (A-1-1). The arrays are filled as follows: 
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SR(1) - SR(X,0) SM(1) = SM(X,0) 
SR(2) = SR(X,1) SM(2) = SM(X,1) 

(A-1-23) 
SR(10) = SR(X,9) SM(10) = SM(X,9) 


ο 900 on I = 1,10 controls the columns. The 


columns of submatrices X, and X, 


ι acion counted by DO LOOP 1000. In each step in 


are filled in sequence for 


LOOP 900, four positons in submatrix X, ere filled with the 


1 
ΘΕΑ. four positions in X, аге zeroizeg. The real 
ШІ ОТ the applicable integral is placed in four positions 
in X. via the dummy variable SUMR. Similarly, the imaginary 
pare 1s placed into X, ΙΙ Сапге 15 (акеп to obtain 
EM ΞΕ--.515. Όποε the LHS is set for each station, 
four rows in the RHS matrix BB are set according to (A-I-22). 
The station is then augmented by the interval, XINT, and 
the procedure is repeated пала positions in AA and BB 
are filled. | 

ПО т СастсиЕ the matrix prior to solution is printed 
G ari EB 15 destroyed in the solution process. 

Προ πα 15 Solved for the unknowns using library 
subroutine DGELG. The inputs are the coefficient matrix, 
RHS matrix, coefficient matrix size (here set to 40 by M=4L), 


RHS matrix size, and a significance measure (here set to 1.0E-6). 


The solution is returned on BB and any error indication on IER. 


25 





ЕР Ш/ о BB(IO) Contain the real parts of the coefficients 
and BB(11) to BB(20) contain the imaginary parts. The parts 
EEcombined in complex form and placed in complex array ZB. 


he form of ZB is 


πλ DE ο ο = 
A O E = 
ZB(11) = DELT(1,0) = (A-I-2H) 
ZB(20) s DELT(1,9) = 


Finally the desired coefficients are printed out on 


Computer paper output and on cards. 
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 A-II. SUBROUTINE WSUM 


The purpose of WSUM is to provide the control for the 


computation of the integral 


X N ο 2 
EEE E е. 
Т = дух D Q τ 
Deets demonstrated in Section III-B, (A-II-1) may be 


reduced to 


ο... s _TT_ 
1. 55 (8 ius uy (А-11-2) 
мһеге 
ς 
n s Е P o exea 
NR 52, G-DTQ-FDT 5 
p (ІІ-30) 
COS [2KAU+E/U) | dU 
and 


М41, QQJ-lQQQ4N-J41 ΧΑΚ - ς 
την W ы 

ja DR | 
CII-31) 


SIN [2KAU+E/U] ) aU 


where che fixed parameters A,B,C,D, and E are set in MAIN 
[Refer to Eqns. (А-І-11, -12, -13, -І1!, and -15).] and 


carried to WSUM in common space HOLD. 


ШІ 





aen the integral values, the summation is a simple 
computer manipulation. The solution of these integrals 
involves three variations of the same basic numerical inte- 
gration technique. The numerical integration of (III-30) 
and (III-31) is accomplished by a "pyramided" Simpson's 
rule using end correction. This technique is presented in 
Appendix A-III where Subroutine ZSUM is described. ZSUM 
Ires the numerical solution of the basic integral for 
СУБУ values of J,K,A, and upper and lower limits, XO 
and XA, are specified. 

The difficulty encountered vith the numerical integration 
of this basic integral occurs as the result of the lower 
ie, AA = 0. For J < 2.5, the amplitude of the integrand 
Mey be unbounded near XA = 0. Further, except in the special 
КО еге K = 0, the SIN and COS functions in the integrand 
Preome unmanageable near zero where the controlling Per 
in their arguments, E/U, becomes large. The resulting 
infinitesimal period cannot be treated by numerical tech- 
niques. For a M As One may neglect the troublesome 
КОШО ОП of the integral near zero. This is permissable 
since the exponential function of (-B/U) forces the ampli- 
tude bs σον ο ος Aintervel where the integrand is ill- 
behaved. When A = 0, for J < 2.5, one may not neglect the 
integral near zero. Γον υπο ο P tiie Case, one may 
approximate the ill-behaved portion of the integral by 


analytical methods. 
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Eu - 0, the basic integrals in Eqns. (III-30) and 


(III-31) become 


ELE 


c 
SUMR = /  UYTéCOS[DU + E/UJAU 
0 
апа 
ΤΙ ΞΞ 
SUMI = ; к + Е/0]а0 
0 


(III-33) 


(III-3/) 


iimroancing the correction term of Section III-B, one may 


write (III-33) and (III-34) as 


τὰς 
COC + 0“ 72 соѕ[ро + E/UJdU 


SUMR f 
AL 


and 


X*l , & 
COS *  [  U'""*SIN[DU * E/UJQU 


AL 


il 


SUMI 


where the correction terms are 


AL 
coc = U 
0 


S 
J-g 


OSEE ИЕ 


and 


AL 
τη 
0 


_ 
J-ESIN[DU * E/U]dU 


ο. 


(A-II-H) 


(A J ο 


Rec 


veras COC and COS are prepared by Subroutine COR. 


The technique as presented in Appendix A-V provides an 
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approximation which reduces the correction integral to the 
well known Fresnel integrals. 

The limit AL is selected such that (1) the numerical 
integration technique may be applied to the second term 
in (A-II-3) and (A-II-l), and (2) the correction approximation 


Nu ud over its interval. 


ΓΕ ANALYSIS 

(Refer to the Subroutine at the end of this Appendix.) 

EuDPoutine WSUM is called to provide the block of L 
integrals at each station on the blades where the boundary 
Eusvlons are enforced. (In this investigation L is set 
at ten.) The subroutine inputs are DK = K, DL = A, RHO = R, 
Em X. 

In the first steps several fixed parameters are estab- 
lished for convenience. EPS = 1.0E-8 is the convergence 
Parameter used in the numerical integration process in 
routine ZSUM. LIM = 50 is a step limiter also used by 
ZSUM [Refer to A-III]. PI = T is used in later calculations 
in WSUM. 

test Is made to determine which of the three cases 
Eeen dealt with. If K = 0, branching is executed to 
Step 1. If not, the next step is encountered. If K # 0 
and À = 0, the routine branches to step 2. The case K # 0 
and À Z 0 follows at step 3. 

EINER = 0 

The intention here is to use the numerical integration 


process over the entire interval. Since no SIN or COS 


104 





Ше iOS enter, one may. approach closely to zero at the 
lower limit (XA is set to 1.0E-8). The upper limit is set 
EN 0 - X + I. Subroutine ZSUM is called to compute the 
ES J - 1,2,...,10. The integral in (II-30) is 
returned on SUMR and the negative of the integral in (II-31) 
pa ctcurned on SUMI. (Note that in this case SUMI = 0.) 

Env slues are stored in the vector arrays AC and AS. The 
routine then jumps to step 1000 where the summation 


gnd@eated in Eqns. (II-30) and (II-31) takes place. 


— 
== 


б гер 2, λ 0, K # 0 


Ша са 5 step the special case of À = 0 is computed. 
The Fresnel integral is obtained, and where applicable, the 
Semrection is computed and summed with that portion of the 
integral obtained by numerical methods. 

Observing that the value of X does not affect the 
correction integrals, they need be computed only once: 


Непсе а test is made: 


Eye Gein) ec Ome Oar 
Step 22 bypasses the correction computation. 
The first steps in computing the correction integrals 
includes the preparation of inputs ο -Ἠτεπττπε СОК and 
Subroutine Red. The rationale for these functions is given 


in A-IV and A-V. They include the following: 
AL = RHO/100 (A-II-7) 


E/AL (А-ІІ-8) 


AR 
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- 


55 


ST CENE) (A-II-9) 


ce COS AD) (Л-ІІ-10) 
Hocohkbpesuel integrals for upper limit AR are obtained 

' Hbrary Subroutine CS. The explanation for CS is 

ЕЕЕ оп call from the library and а copy is found at the 

end of this Appendix. The only argument required is AR. 


INéStollowing is returned. 


AR 
cec -| 505) ат πω. 


үгтт 
A 


R 
sss =| SINC) an ΠΤΙ) 


E vant 7 


Mampmecparation for their use in Subroutines COR and RED, 


(A-II-11a,b) are modified as follows: 


σος - = [5 - ccc] Ее. 
апа 
SSS = ы. [5 = sce : (Оор) 


Hefer to A-V for an explanation of these functions. 
The various inputs in hand, Subroutine COR is called to 
Pocpzrersche block of L - I0 correction integrals for Eqns. 


E55). The arguments of COR are AL, D, E, CC, SS, 
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COG; апа 555. Тһе value of (A-II-5,-6) for each J is 
femurned on CCOR and SCOR respectively and stored in arrays 
coc and COS. COC and COS are of vector size ten and are held 
in common with MAIN for permanency. Also note that (-SCOR) 
' τοτεὰ in COR for compatibility with the numerical solu- 
tion of the remaining integral. The correction terms are 
printed for reference. 

Finally the remainder of the integrals are computed 
numerically in the same manner as was used in step 1 but 
with lower limit AL. The result is summed with its correction 
BEusMScored in AC and AS. 

Step 3. K # 0, A # 0 

In this step, one relies on the exponential function to 
Meauce the amplitude of the integrand near zero. That portion 
спе integral from zero to 1.0E-3 is neglected. The 
ως. cal solution is called as before and the results 
Seemed in AC and AS. 

ШЕР сте Values of the integrals for J = 1,2,...,10 are 
in hand, the indicated summations, Eqns. (III-30,-31), 
must be carried out. Step LOCO ερ ΤΠ орегайтоп. The 
1ο} is carried out іп а straight forward manner for 
апа І 


ЕП ОТ the L = 10 integrals. I are computed on 


NR NI 
SUNG and SUMS respectively. (Note that IFAC(M) is a function 
subroutine generating the factorial of M. Refer to A-VI.) 
Meco. leteda integrals are stored in arrays SR and SM and 


r j I ed Tor reference. 
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Recalling that equation (A-II-2) contained Q = A+2iK, 


a final complex manipulation must be made. 


where 


Then 


where 


and 


for 


Finally 


I 


3 


Б 


Q = 


о = 


JOE oS ο... 


ARCTAN(2K/2A). 


RIQ _ : 
Бү = ZÍBK + 1AK] 


AK 


BK 


. 7 
,R [η 


2 ym 


SIN(G) 


COS(G) 


NR 


ARCTAN(2K/A). 


(A-II-2) may be written 


= = 1 T 
2 (IypBK TAX) + 12( Ты БАК. d 
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NI 


Observe that 


Fl). 


(А-ІІ-13) 


(A-II-14) 


Cane 


(o 


(A-I-18) 


(А-І-16) 


(A-II-15) 





ιν end BK are compuved in MAIN to eliminate 
τση and are carried to WSUM in common. DO LOOP 1030 
computes the real and imaginary parts of (A-II-15) for the 
ten integrals and returns these values to MAIN in arrays 


SR and SM. 
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A-ITIT. SUBROUTINE ZSUM 


subroutine ZSUM computes the value of the integral in 
the general Eqns. (III-30) and (III-31) by numerical 
methods. 
XO 


ν᾿ 
ХА 


J 


E 
SUMR 2 PPI CU EU ICOSIDUFE/U Т ОШ М (ATI 


XO 
о 
ХА 


J 


_ 
SUMI ZEXPLCU-B/U JSIN[DU+E/U JAU  (A-III-2) 


HnEEULechnique utilized is Simpson's rule with end correction. 
Noting the behavior of the integrand, especially near 
ШЕ 02 15 185 not satisfactory to compute the integral using 

a fixed interval. Instead a "pyramided" scheme is used 
which automatically reduces the size of the subintervals 
until the scheme can compute the integral to within the 

watSs ОТ а Satisfactory convergence criterion. 


i nceceneral form of Simpson's rule with end correction 


is 
EE I 
SUM = 15(7Ү,%16Ү%18Ү,%16Ү,%...%16Ү, 147Ү,) 
2 
H ! : 
+ i5 14 -YN ) 
or -1 N- 
H 2 2 y2 
ze + RR z 
SUM i5 C1 (94 * Y, 16 Z Y 112 Yous) trey Үн ) 
і-1 1-1 
(A-III-3) 
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t 
where - and X; Pe πας Value vot place integrand and the 
derivative of the integrand evaluated at station ad 


Hees the length of the subinterval [Refer to Fig. 8]. 





Let Xie2 


Figure 8. INTEGRATION DIAGRAM 


~ 


Heure 9 demonstrates the "pyramid" algorithm. X4 and 


Xx sache lirst and last station of the interval over which 


Y(X) is to be integrated, i.e., X4 = XA and Xy = X0. N is 
the number of stations at which the integrand must be 
evaluated. There are (N-1) equal sub-intervals of length H 


ШЕТІ = Xy - X4. 


Step N H Xy сое j [O οἱ As EN 





ΙΤ |INT/16 |7 16 14 16 14 16 14 16 14 16 14 16 14 16 14 16 


PYRAMID SCHEME 
ΕΙ Τη 
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Euecaecuracy of Simpson's approximation inereases with 
decreased H. Тһе pyramided technique progresslvely computes 
the integral, doubling the number of' sub-intervals each time, 
until the value of the computation ceases to change to 
Es some desired criterion, EPS. 

Figure 9 shows the change in N and H with each step. 

' Defficients for each step are shown in their relative 
Fo l ions according to station. Observe that once a Y: ШЕ 
ΙΙ 64 for a given step, it is used in each succeeding 
step and need not be reevaluated. Further note that the 
first time a Y, is used, its coefficient is 16; thereafter 
Іс 11, 

The integrand in (A-III-1) is 


2o gl 
4549 m U 2 


EXP(CU-B/U)COS(DU+E/U), (A-III-!) 


and it may be shown that 


! | 7 
Y (U) =U? -ZExP(CU-B/U) | COS (DU+E/U) (J-5/2+CU4B/U) 


R 
(A-III-5) 
-(DU-E/U)SIN(DU+E/U) } 
Шие тпсергапа іп (A-III-2) is 
5 
Y (U) = UJ "ZEXP(CU-B/U)SIN(DU*E/U) mm ο 


апа 


"πι 





- 


t -7 
Yr eu = 07 FEXP (CU-B/U) { SIN(DU+E/U) (J-5/2+CU+B/U) 


(A-III-7) 
+ (DU-E/U) COS (DU+E/U)| | 


Seer ANALYSIS 
(Refer to the subroutine at the end of this section. ) 


Subroutine ZSUM arguments are as follows: 


XA = Lower limit 
XO = Upper limit 
J = Summation index set in WSUM.[Refer to Eqns. 


(11-30,-31).1 

The following parameters are brought to ZSUM via common 
ШОО ШОГО: D,C,D,E,EPS, and LIM. B,C,D, and E are computed 
in MAIN [Refer to Eqns. (A-I-12,-13,-14, and -15).]. EPS 
and LIM are set in WSUM. XA, XO, and J are also set in WSUM. 
The values of (A-III-1, -2) are returned on SUMR and SUMI. 

Equations (A-III-4,-5,-6, and -7) are first evaluated 
at the left endpoint. AA = CU-B/U is the ent of the 
exponential function. If AA is less than (-172), it is 
sn d char the exponential term drives the integrand and 
its derivative to zero by branching to statement 6 (Computer 
underflow errors are generated for smaller values of the 
exponent.). BB = DU+tE/U, the argument of the SIN and COS 
functions, and CC = DU-E/U, a term used in Тү; are computed. 
The evaluations of Y p (XA) and Y (XA) are obtained on YAR 


and YAI from the arithmetic statement functions YR and YI 
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where YR is (A-III-!) and YI is (A-III-6). EE ana FF are 
EuUEonal terms used in the derivative. Finally the values 


1 1 
of Yn (XA) and Y, (XA) are computed on YARP and YAIP. 


I 

Statement 6 sets all values Y (XA), Y (XA), ου» 
апа Y, (XA) to zero for small values of the exponential 
function (those cases where AA < -172). 

esame process is executed for the right end point. 
Here Y,(XO) - YOR, Y.(XO) - YOI, Y, (ХО) = YORP, and 
Y. (XO) - YOIP. 

[1555 step in computing the internal points is to 
m interval of integration in DINT = XO - XA. Next, 
S Summation parameters are zeroed including SOR апа 
See che sums of the Үр, 5 апа Үт: 
by the coefficient 14 [Refer to Eqn. (A-III-3)] and SER 


's which will be multiplied 


and SEI, those sums to be multiplied by 16. SP sets the 
number of subintervals in each step (initially set at 2) 

and K sets the number of terms in SER and SEI to be computed 
each step (initially set at 1). 

DO LOOP 100 controls the "pyramid" level or step from 
one to LIM. His the length of each subinterval in each 
ЕВН = DINT/SP. 

Ón each step, the Y,'s computed the previous step for 
coefficient 16 (summed in SER and SEI) are added to those 
with coefficient 1! (SOR and SOI). SER and SEI are then 
Zeroed in preparation for the next computation. XE is the 
E nm vardiable in this portion of the calculation. The 


initial Y; is computed at XE = XA + H. 
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DO LOOP 10 computes the K new Y,'s τον; επα SEN. 
AA and BB are computed as before. Yp (XE) and Y (XE) are 
Euputed on YEH and YEI and summed on SER and SEI respectively. 
Finally the dummy variable is augmented by 2H and the loop 
Semoinues. 

SUMR and SUMI hold the values of (A-III-1) and (A-III-2) 
Гог any given step. TEMPR and TEMPI are used as temporary 
storage for the values in SUMR and SUMI computed in the 
previous step. Once a new set of values are computed in 
БЕП апа SUMI, they are compared with the previous values. 
"EN rh do not differ from the previous values by more than 
eee onvergence criteria is fee and the loop terminates. 
сарша сле process reach LIM iterations without convergence, 


a warning is printed and the existing values of the integrals 


are accepted. 
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A-IV.  SUBROUTINE COR 


Kur Line COR controls the computation of the correc- 
Eu usintesgrrals (A-II-5) and (A-II-6) for the case À = O. 


The integrals are 


AL J-£ 
CCOR = f U 4COS(DUTE/U)QU (A-IV-1) 
0 
апа 
АТ, 1-5 
SCOR = I U 2 SIN(DUIE/U)dU (A-IV--2) 
0 


- 


where D and E are specified in (A-1-14,-15). An analytic 
Exmation for these integrals exists in theory; however, 
Bempurer numerical error prohibits reliable use of this 
πε for J > 5 and R > 2 [Refer to the Discussion of 
Results for an explanation]. For this reason, the approxi- 
mation is used only for J < 4. This range includes those 
troublesome sss having integrands with infinite ampli- 
tude near Όρο πήρα Stangard numerical integration is 

used for J > 4 where the amplitude of the integrands becomes 
small near zero thus enabling one to approach Zero sufficiently 


οσο το ensure that the omitted interval is negligable. 


CASE FOR J < Ц 
Under the assumption that, for a properly selected AL, 


as U approaches zero the term E/U dominates DU in the 
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arguments of SIN and COS then equations (III-40) and 
(I11-41) are reasonable approximations of (A-IV-1) and 
(A-IV-2) even when a reduced number of terms is included. 
ENS selection of AL is AL - R/100 which ensures 
(E/U)/DU > 1.0E4 for U « AL. 


Subroutine COR carries out the operation of (III-40) 


and (11-41) where 


AL +, € AL y 3 
CCOR = f U'"£cos(E/U)dU - D f UP" ZSIN(E/U)dU 
0 0 
ο... СШ 
-5 у 0972соѕ(в/0)ау + = f wW* sınE/u)au 
0 0 
(A-IV-3) 
and 
AL Ίο Tu Mm 
SCOR = fs UJ £ZSsIN(E/U)dU 4 D / U/-?cos(E/U)qaU 
0 0 
О сү en 
EMI S UN dU a cos (E/U) au. 
° ϱ 6% 
(A-IV-l) 


шего сев пассигасу та спе computation of the integrals 
ПОШ I TOL IO RED, a power of U greater than 1.5 is avoided 
by reducing the number of terms in the truncated series by 
one each time the series is used. That I for J = 1, Tour 


terms are used, for J = 2, three terms are used; and 5ο 


forth until only one term is used for J = H. 


1ο. 





ioe FOR J > 4 | Ξ 

ШӘП thas case Subroutine ZSUM is used to numerically 
116 the integral over the interval XA to AL where XA 
NN "cked as close to the origin as computer time and 
 Шрасу warrant (XA = 1.0E-H was seclected for these 
results). It should be noted, however, that in some cases 
maemvaiue of the integral over this interval may be less 
ШЕШ (Пе convergence limit, EPS, used in ZSUM in which case 
Paiemtcosults may suffer from inaccuracy. Provided the error 
is on the order of EPS (EPS should be the upper limit) such 
E cuPacies are not significant and should not affect the 


end result. 


eee ANALYSIS 
(Refer to the subroutine at the end of this Appendix.) 


ime anputs to Subroutine COR include: 


ШІ Upper Limit of thewecorrection integral 
J = Summation index in equation (III-30) and (III-31) 
2x? 
um 2-2 
AS + HK 
E _ KR 
2 
CI COS (E/AL) 
S = SIN(E/AL) 
AL, 
SS = / —, SIN(E/U)QU 
0 ША 
AL 4 
cc = / — COS(E/U)aU 
y ue 


« 


126 





D and E are fixed parameters computed in MAIN; C, S, 
55, and CC are formed in Subroutine WSUM. CCOR and SCOR 
return the values of ο and (A-IV-4) to WSUM. 
meee thst step is to determine if J= 1,2 « J « h, 
or J > 4. For the case J = 1, four integrals are computed 
in τς or the form 
AL 


CINT = f U 
0 


Mei | 
2 COS(E/U)AU η ΤΟ, 


апа 


AL М-5 
τς U STNCEAU)QU (A-IV-6) 
0 


Еке М = J, 741, Jt2, and J+3. These integrals are easily 


1η four step DO LOOP 10 utilizing Subroutine RED. 


RED reduces (A-IV-5,-6) to the general form 


AL 
Integral - f(AL,N,E) & 7 on 


0 US 





(E/U)aU. (Beyer) 


The values of these integrals are stored in two arrays, 
СОКЕ cot damension four; Note that should a full four 
term truncated series be desired for J > 1, only one addi- 
tional integral is needed in (A-IV-3,-4) each time. That 
terms ot the J of the preceding calculation, values 
ο ΝΕ ana SINT are required with M = 741, J+2, J+3, and 
J+4, Since the first three are already on hand, DO LOOP 


EN still) into SI(TI) aná CI(1+1) into СІ(І). Тһе 


Ar 





fourth position is filled with the newly calculated values 
of CINT and SINT for M = J+4. However, assuming error in 
ШРШ Тос J > 4, it was elected to truncate the series by 
one additional term for each case, J = 2,3, and 4. Hence, 
Eus set into 51(4) and CI(!). 

Encossumes that this reduction in accuracy is still 
Eccprcable compared to the error introduced by an attempt 
IE Ue numerical integration. 

In the case J > 4, a branch is provided to statement 21 
iets subroutine ZSUM is used to calculate the integral 
between AL and XA. Although this procedure simplifies 
the handling of COR results by WSUM, inaccuracies of the 
order of EPS are expected as previously mentioned. One 
Шаа ехресс greater accuracy using another technique 


"Ου ο not compute the integral piecemeal. 
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SUBROUTINE СОҚ (АС,О,Е,4.2С,5.55.СС,5СОҚ,ССОК) 
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ZSUM(CXA ,X0, J4 CCOR, SCOR) 
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A-V.  SUBROUTINE RED 


oubroutine RED provides the value of the general integrals 


А J-E 
SINT ES NAO AU (A-V-1) 
0 
апа 
А J-£ 
Cie А ше AU QU (A-V-2) 
0 


КОО ОС пр them to functions of the form 


Ay SIN 
SINT = F,(A,J,E) + e = Sag (E/U) AU (A-V-3) 
and 
А SIN 
ο πο ο ο r (E/U)QU. (д-У-!) 
2 211% С05 


The integrals in (A-V-3,-!) are funetions of the Fresnel 


inteprals as follows: 


EB E ο {οποι 
сс = у 2 cOS(H/U)aU = [E [I - C(E/A)] (A-V-5) 
0 9% 
апа 
А 1 E gar 
SS = f —-SIN(E/U)dU s J—[z - SCE/A)] ne ΞΕ, 


lol 


where _ 


E/A 1 
/ 
0 V2TT 


ӘСЕЛ” 


CoS Gyan Каат) 





апа 


Е/А 
| T 


S(FZA) 
0 V2rT' 


ΙΙ Dal: (A-V-8) 





That is, SINT and СІМТ for any J may be reduced by 
Ширевгагіоп by parts to functions of the basic Fresnel 
ΙΓ ΕΡΡΟΙΘ. 


Integration by parts provides the following relations 


for (A-V-1). 
ШОО = 1, 
А 
NU e COS(E/U)dU. MA=V=9) 
OS 
mor J = 2, 
A 
CINT = E 


2[A*COS(E/A) - E f -SIN(E/U)AU]. (A-V-10) 
000 


шісі тог J > 3, 
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CINT - CEPR [α΄ . 


- (22-30! 
J o. 
[z4.1] 
° ΄ N-1,4N- ON Е 
-[COS(E/A) x Е ον ο Am 
N=1 Scout 
ce 
E NO NR 
- SIN(E/A) Y [2(J-2N-1)]1(-1) 7-129 7N-1,J-2N-1,2N-1 
TESI ENZO 
y тит} (Дегені 
where 
A 
"RU 
BENEI a) 
or 


ο... 
Int = 220-3;7-1( 1)6%); =, SIN(E/U) GU, AO 
0 U 


(A-V-12b) 


Similarly, it may be shown that for (A-V-2) one obtains: 


Fon J = l; 
P 
SINT = / L SIN(E/U)dU, (A-V-13) 
DEUS | 
ifJ-2, 
1 A 
SINT = 2[A^SIN(E/A) * E / +, 008(E/U) au), (A-V-14) 
0 U 
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гот J > 3, 


το Vy. 
Ju 
[4.1 | 
| [sinc ED C D ын) 
| N=1 MES 


j 
L. 
2 N-1,N-1,J-2N-1,2N-1 
"τ; -- A 2 | 
ME PEOPLE E 


(A-V-15) 
+ тнт} | 
where 
ην [5] Ay J 
INT = 2 Е 730-1) < / Jy BIN(E/U)dU, (-1) «0 
000 
η a 
Or 
j 
(A 
Шз (y σα m 0. 
ο A 


CES Sle) 


moe notation LARG] indicates the sign of ARG with the 
largest integer less than the absolute value of ARG. The 
use of [ARG + 0.1] ensures the proper value is selected and 


Bundes the possibility of computer round off error. 


ӘШІР ANALYSIS 
(Refer to the subroutine at the end of this Appendix.) 
The arguments of Subroutine RED include: 


ТЕРЕ (тыс ol the Correction integral. 
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E = KR“ /2 


- 


Je =) Related to the summation index in (III-30;-31) 
Scc SIN(BZA) 
ПЕЕ Сос СЕХА) 
а 
OD = И y SIN(E/U)aU 
2 
O U 
EE 
СС = f —, COS(E/U)AU 
o y 


Bis set in MAIN and carried in common to WSUM and as 
me Όπου through COR to RED. A,S,C,SS, and CC are 
I ir їп Subroutine WSUM and are carried as arguments to 
RED via COR. J is set by COR. Тһе values of (A-V-1,-2) 
Es ccurned to Subroutine COR on SINT and CINT. 

EM coc 1S made to determine if J « 3. If it is, the 
une branches to statement 1000. 

If J > 3, equations (A-V-11) and (A-V-15) apply. Τι 
a o are the summation limits where INT(ARG) = [ARG] 
deseribed previously. The primary series DO LOOP computes 
the first summation in (A-V-11, -15). The secondary series 
DO LOOP eomputes the second summation. А test is made to 
check the criteria of (A-V-12,-16) and proper branching is 
initiated. Statement 31 computes (A-V-12a,-16a) for 
e < 0, and statement 32 computes (A-V-12b,-16b) for 
eL) EET dy. a statement ЛО combines the terms to 


form (A-V-11,-15) on SINT and CINT. 
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Ш пе case ς ο, -ᾱ test is made to further restrict 
the case to J = 2 or = 1. For Ј = 2, (A-V-10,-14) are 


used. For J = 1, (A-V-9,-13) are computed. 
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SUBROUTINE ВЕО(А, Е. Ј, 5, С. 55, СС, 5 ІМТ, СІМТ) 
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AVI, SUBROUTINE IFAC 


Function subroutine IFAC returns the value of the 
argument factorialed, e.g., IFAC(5) = 5-4-3-2+1. It has 
provisions for IFAC(0) = 1. A copy of the subroutine 


оон. 
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APPENDIX B 


ЭЛ ee POENTIAL SOLUTION 


B-I. COLLOCATION TECHNIQUE 


in order to investigate the perturbation velocity 
potential along a cascade blade for the steady state and 


m phase shift in adjacent blades, one may use the following: 





X X 2 
1 1 i АВ 
D ος E = 5% f ——XP[- κε] 
pes VEMM CD υπ I 


X aN — 
ШЕ. S AR 
crop ji EXP[- Jas 
Í: πετ = 
Vix Eon EIE m 4025) 





X QN 


f аз (III-lT) 
О 





The second and third integrals are of the form 





X aN 2 
T = f КОР BA E аз, where B= AR 7A. 
ps | (ED 
Letting 
Us OS» (B-I-2) 
one may write (B-I-1) as 
κ ote а 7 ο το 
3⁄2 
B U 
xt] 
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OBS 


] i M 
ШИ A ERU απ. (B-I-4) 
ы 
В U 
ЖЕРДІ 
The binomial expansion provides 
Gesu)’ = 3 0-10) м (виху)К, ee) 
ке (N-K) IK! 
Hence (B-I-H) becomes 
ET со 7 
2 δα, (-1) B -(K*z) _-U 
| ШЕРІ V ang NC EREXE Б "ERU e GM (B-I-6) 
AFL 


It may be shown that integration by parts provides 


ki -X _ T e 
πι ο ο ο πι πετ ο 





n 
п 


ο Eo EO Cela es) 


1 
M-1 


where M may be any half odd integer greater than zero, 
emer = 1/2. 3/2, 5/2, etc. Continuation of the integration 


Process yields a general relationship 


М-1)! 


ME -pé n-5 1-1). 
cr Ë Σ e 
p2 X LLAMA 
1 
I MO AT ERE (E) (BESOS) 
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where 


il 


U 
„5 &- dü = —[1 - erf(F)]. eo, 


т үт 


ERFC(F) 


li 


For the case N 0, (B-I-6) is 


Е 
lI 


VP / EU 
z dU. (B-I-10) 
xau 





The last integral in (III-47) may be written 


X oN N „N-K 
! (X+ 
M SE кп. E wie E ΝΕ e: n Garn) 
E = S 


and observe that for the case N = 0, (B-I-11) becomes 





Е s? (o 
1 


Combining Eqns. (B-I-6,-10,-11, and -12) in (III-47) 
aud replacing N with J-1 and K with K-1, one arrives at 


the final computer compatible equation 


II 


(26-199 xn + 2/8 SSUM[-3/2, b/(X+1)] 


ТА 
E Фәса лмаг, ЖО түш 
., κ υπνο ROBES —| 0,2-1 2K-1 


K-1 1 
=D, 5.18 Y8 SSUME-(Kt5) ,B/(X*1)1] 1 (B-I-13) 


ko 


> 
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where 


SSUM(A,FSQ) = 5 c Ge 
FSQ 
ENS computed in Subroutine SSUM and will be dealt with 


in Appendix B-I-B. 


B-I-A. THE MAIN PROGRAM 

Ihe purpose of MAIN is the computation of (B-I-13). 
Inputs to the program are the controlling parameters, blade 
spacing, R, and thickness parameter, X, as well as the 
n S OI the interference coefficients obtained from the 


eommeocation coefficient solution of Appendix A. 


\ 


Step Analysis 
(Refer to the routine at the end of Appendix B-I.) 
Double precision is used in Subroutine SSUM and single 

precision is used in MAIN. Double precision arguments to 

and from SSUM are carried on variables beginning with W. 
ШЕ ο μυς πλ А, ΡΙΞΕΙ απ ЕНО - KR, and the 

coefficients of the upper and lower interference functions 

(0 апа 7 


lm ІШПЕ 
data cards and printed out immediately for verification. 


m= 1,2,...,10) are read into the routine from 


(Note that K, the reduced frequency, is used for accounting 
purposes on the output only.) Array DEL(20) carries the 
coefficients with the lower blade coefficients in the first 


Leo positions. 
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s puts positioned, various fixed parameters are 


περα and the station, X, is initialized. 


B = ЕРІ 
BH = JB | 
С = 1/утх 
Dos уо 
X = -0.9 


DO LOOP 100 is the station control loop which generates 
EN crential value at each of 20 positions on the blade 
wee -0.9, -0.8,..., 1.0). 

Іп preparation for the use of Subroutine SSUM, the lower 


ΠΠ 15 formed, 
WFSQ = B/(X+1). 


КОШОМ 15 called to form 


о -U 


SSUM[D,B/(X+1)] = $ | ez dU. 
U 


—_ eae 


B 
Al 


With AA = yX+1, the first two terms in (B-I-13) are formed: 
BB = 2/X+T + 2/B' SSUM[-3/2,B/(X+1)J]. 


DO LOOP 200 accomplishes the J summation in (B-I-13) 


while LOOP 300 carries out the K summation. (Recall that 


1545 





το.) Each term in the summation is built up in a 
Seraitghtforward manner. Subroutine SSUM is required to 


compute 


paro a 


— 


SSUM(A,WESQ) = > 1980 


where 
- 1 


mume t lon subroutine teat coma SOmused eels) subroutine 
55 the factorials of its argument and is explained 
in Appendix A-VI. 

ive summation loop completed, the summation term is 
combined with BB and C to form POT - Vo (5,0) of equation 


(B-I-13). Finally POT is printed together with its station. 


% 


B-I-B.  SUBROUTINE SSUM 
ime purpose of SsUM is to produce the value of the 


general integral 


со 


аа aa (Β-1-15) 


- 


nro| 


p? 
Only the derivation of the procedure for A < O will be 
presented in detail here. The case for A > 0 is derived in 


similar: manner. This Matter case is not of interest in 


ше potential solution. 
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lucsoneline funetion subroutine, INT(A), generates the 
sign of A with the lareest перо ο επ μας СЕУ 
INT(-3/2) = -1. The subroutine, ERFC(F), produces the 
ement of the error function using pS as the lower limit. 


It should be understood that (N/2)! - > (5 - 1): (5-2)... 5 


and that (-1/2)! = 1. 
Now for demonstration purposes, (B-I-8) is expanded 


lx e cases A = -3/2 and A = -5/2. 


A= - 3/2 
т s^ Πεμ» 
Bf 255 4 = > = -2УПЕНРС(Р)1 πο, 
F е 
ШЕ 5/2 
1o ge^axk-l gpl 2... a | 
> JE > s MEDIE = 
Ex „Fr 
(B-I-17) 


J routine may be generated in the following manner. 
Observe that if one defines I = INT(A), there are I terms 
in the series (Those terms not e the error function 
term.). The numerator is increased by two іп each succeeding 
term. The denominator of the first term includes (21-1); 
each succeeding denominator in the series is multiplied in 
turn by (2I-1)-2, (2I-1)-4, .... The power of F in the 
denominator begins with (2I-1) and is decreased by two in 


each succeeding term. The sign changes with each term. Also 
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note that the numerical coefficient of the error function 
ШЕГІП 15 the same as that for the last series term with the 


sign reversed. 


Step Analysis 
(Refer to the subroutine at the end of Appendix B-I.) 
Ihe inputs to the subroutine are A and FSQ - pS. The 
meme. or (B-I-15) is Fe Liy q on UM ІШ the first steps 
various fixed parameters are formed, the summation parameter 
is zeroed, and INT(A) is determined. A test is made to 
Ec Пе 1Г А > 0 or A < 0. Here only the branch A < 0 
memoreincerest although the structure of the unused portion 
eerie subroutine for A > O is similar. 


ШЕШІ < 0, statement 10 tests for the special case 


A=5 (or I = 0). This case cannot be generated by the 


algorithm used for other cases. If A = - 5 , Statements 12 
and 14 set 
SUM = LE EREC(P) 


and the routine returns to MAIN. 


l 
If A 752 


several preliminary steps, DO LOOP 11 forms the series 


statement 13 begins the algorithm. After 


pointed out previously. The method takes advantage of the 
Doservacions made in the introductory material. The series 
mermed on SUM, it is divided by the exponential and the 
Une cion term Is added. 


The result is returned to MAIN. 
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Ert. 


Development of the Sandeman transonic 
interference solution for the steady case 


showing perturbation potential function 





X 
w(X,0) = - d. к ш as 
V TÀ | 4 У Х-5 
со ye 
+0 VS) ΕχΡΓ-λ 
A Е V X-S 


compare this solution with the cascade 


Exod to blade spacing R by d = 


Further, as shown in Egn. (III-18), 


v(S) = Y (X,Y = 0) = -l, for K = 


SONIC WIND TUNNEL WALL INTERFERENCE TECHNIQUE 


wind tunnel wall. 
resulted in the 


ОПСС оГ = 


ма. T 
Х-5 
(III-60) 
problem, (d) is 


R/2 [Refer to Fig. 7]. 


D. 


EONollows that the potential may be written: 


i г 1 
οσο) fe --- E MIS 
Y TA -l v X-S 
X 272 
a0) —— EXP[- А es] ds 
в 


Iinteerating the first term and making the 


in the second term one may write 


со 


со 22 
1 AN^R 
%(Х,0) - --- {2Υχ11 +2) 
VTA N=1 on 
ΠΟΤ 


>> 


S 


(B-II-1) 


proper substitution 


-U 


(B-l) 





The basic integral within the summation mc cc τη 


Mois Oi the error function, i.e., 


oo E р 
е U 2 i 
К dU - ——— — — 2s ERFOCR) MEE 
αυ F 
F 
where 
222 
2 XNER | 
Р = тС) .: 
апа 
ВОВСЕ) - 3 - erf(P). (B-II-5) 


ος the potential function: 


ψίχ,0) Л т ὁ ὧν ΠΠ 
ш N=1 


АМ 


5—5 | 
?R^  gRFC[ T | (B-II-6) 


-TAN 
Computation of this function is easily carried out by 
τες, The series is truncated when the inclusion of 
additional terms only changes the sum by a value less than 


„Ше convergence criterion, EPS. 


Step Analysis 
(Refer to the routine at the end of this appendix.) 
The inputs to the routine are DLAM = A and RHO = R. 
ШЕ КО 15 inputed for accouting purposes in the output only. 


The output includes the above parameters and the value of 
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НЕ регопграсіоп velocity potential at each of twenty 
stations along the dada 

Л ЖОО ТО parameters and various fixed parameters are 
 Шоптарес prior to entering the control loops. DO LOOP 100 
Кие qes control to calculate the potential at each of the 
Stations beginning at X = -0.9 and increasing by increments 
of AX = 0.2. DO LOOP 50 computes the series summing terms 
 сопуегрепсе і5 reached. Again, the function subrou- 
IEESESEREFC is a standard on line library routine performing 
the operation in (B-II-5). 

Mice the series is computed, statement 51 forms the 
potential on POT which is subsequently printed along with 


wee applicable station location. 
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